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Stokes wave is a finite amplitude periodic gravity wave propagating with constant ve¬ 
locity in inviscid fluid. Complex analytical structure of Stokes wave is analyzed using 
a conformal mapping of a free fluid surface of Stokes wave into the real line with fluid 
domain mapped into the lower complex half-plane. There is one square root branch point 
per spatial period of Stokes located in the upper complex half-plane at the distance Vc 
from the real axis. The increase of Stokes wave height results in approaching Vc to zero 
with the limiting Stokes wave formation at = 0. The limiting Stokes wave has 2/3 
power law singularity forming 27 t/ 3 radians angle on the crest which is qualitatively dif¬ 
ferent from the square root singularity valid for arbitrary small but nonzero Vc making 
the limit of zero Vc highly nontrivial. That limit is addressed by crossing a branch cut of 
a square root into the second and subsequently higher sheets of Riemann surface to find 
coupled square root singularities at the distances ±Vc from the real axis at each sheet. 
The number of sheets is infinite and the analytical continuation of Stokes wave into all 
these sheets is found together with the series expansion in half-integer powers at singular 
points within each sheet. It is conjectured that non-limiting Stokes wave at the leading 
order consists of the infinite number of nested square root singularities which also implies 
the existence in the third and higher sheets of the additional square root singularities 
away from the real and imaginary axes. These nested square roots form 2/3 power law 
singularity of the limiting Stokes wave as Vc vanishes. 

Key words: Surface gravity waves; Stokes wave; Complex singularities of two-dimensional 
fluid flows; Free surface flows 


1. Introduction 

In Part I ( Dyachenko et aZ!||20I6 ), we obtained Stokes wave solution numerically with 
high precision and analyzed that solution using Fade approximation. We showed a con¬ 
vergence of Fade approximation of Stokes wave to a single branch cut per spatial period 
in the upper complex half plane C"*" of the axillary complex variable w. In this paper we 
formulate the nonlinear integral equation for the jump of Stokes wave at the branch cut 
in the physical (first) sheet of Riemann surface. We show that the Riemann surface of 
Stokes has infinite number of sheets as sketched in Figure and study the structure of 
singularities in these sheets. 
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Figure 1. A s chematic of Riemann surface sheets for non-limiting Stokes wave in the complex 
variable ^ (1.6 1 near the origin. The hrst (physical) sheet has a single square root singularity 
at (■ = ixc in the upper complex half-plane C"*" with the the lower complex half-plane C“ 
corresponding to the domain occupied by the fluid. Other (non-physical) sheets have square root 
singularities at <^ = ±iXc. Dashed lines show branch cuts. In addition there are the singularities 
at C = ±i in all sheets which corresponds to w = oo. As well as starting from the third sheet 
there are square root singularities away from both real and imaginary axes at the distances more 
that several times exceeding Xc, i.e. well beyond the disks of convergence |C ± iXc\ < S^c- 


Stokes wave is the fully nonlinear periodic gravity wave propagating with the constant 
velocity c ( Stokes|[l847 1880a I. It corresponds to two-dimensional potential flow of an 
ideal incompressible fluid with free surface. Following Part I (Dyachenko et al. 2016), 
we use scaled units at which c = 1 for the linear gravity waves and the spatial period 
is A = 27T. Nonlinearity of Stokes wave increases with the increase of H/X, where H is 
the Stokes wave height which is defined as the vertical distance from the crest to the 
trough of Stokes wave. Stokes wave has c > 1 and the limit iT —>■ 0, c —?► 1 corresponds 
to the linear gravity wave. The Stokes wave of the greatest height H = H^ax (also called 
by the limiting Stokes wave) has the singularity in the form of the sharp angle of 27 t/ 3 
radians on the crest (Stokes 18806). We assume that singularity of the limiting Stokes 
wave touches the fluid surface at w = 0 and corresponds the following expansion 


z(w) = i-i 

^ ’ 2 




2 ) 


P h.o.t. 


( 1 . 1 ) 


which ensures the sharp angle of 27 t/ 3 radians on the crest. Equation recovers the 
result of Stokes (18806). Here h.o.t. means higher order terms which approaches 0 faster 
than w'^'^ as w ^ 0. Also 


z{w) = x{w) + iy{w) (1-2) 

is the conformal transformation which maps a half-strip — tt ^ u ^ tt, — oo < n ^ 0 of 
the conformal variable 


w = u + w 


(1.3) 


into a fluid domain of inflnite depth —oo < y y yix), —tt y x y tt oi the complex 
plane z (see Figure 1 of Part I ( Dyachenko et a/(]|2016 )). Here x and y are the horizontal 
and vertical physical coordinates, respectively, y = rj{x) is the surface elevation in the 
reference frame moving with the speed c. As discussed in details in Part I, choosing 

z{w)=w + z{w), (1-4) 

with x{w) = u + x{w) and y{w) = u -|- y{w), ensures that z{w) is 27r-periodic function 

z{w + 2tt) = z{w), i (±7r) = 0. (1-5) 
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It was found by Grant (1973) that the corner singularity O) might not be a simple 
algebraic branch point because next order term in the expansion (0) might be a power 
of the transcendental number. Rigorous results on the asymptotics near the crest of the 
limiting wave were found in Refs. Amick & Fraenkel (19871; McLeod|(|1987 |. These results 
were used in Refs. Fraenkel (2007 20101; Fraenkel & Harwin (2010) to construct the 


exact bounds on the limiting Stokes wave and prove the local uniqueness using Banach’s 


contraction mapping principle. More exact bounds were provided in Ref. Tanveer (2013). 


However, the question if log terms in the asymptotic expansion are possible in addition 


to the transcendental power asked in Ref. Amick & Fraenkel (19871 remains open. The 


existence of limiting Stokes wave with the jump of the slope at the crest in 27 t/ 3 radians 
was independently proven by Plotnikov (1982) and Amick et al. (1982). 


In this paper we focus on analyzing singularities of near-limiting Stokes wave. [Grant] 
(1973) showed that assuming that singularity is a power law branch point, then that 


singularity has to have a square root form to the leading order. Tanveer (1991) provided 


much stronger result proving that the only possible singularity in the finite complex 


upper half-plane is of square root type. Ref. Plotnikov & Toland (2002) discusses the 


existence of a unique square root singularity above crests. The existence of only one 
square root singularity per period in a finite physical complex plane was also confirmed 


in Ref. (Dyachenko et al. 2013a) and Part I (Dyachenko et al. 2016) by analyzing the 


numerical solution for Stokes wave. 

We now consider an additional conformal transformation between the complex plane 
w = u + iv and the complex plane for the new variable 


C = tan(|) 


( 1 . 6 ) 


which maps the strip —tt < Re{w) < n into the complex C, plane. In particular, the line 
segment — tt < w < tt of the real line w = u maps into the entire real line (—00,00) in 
the complex ^-plane as shown in Figure 5 of Part I (Dyachenko et al. 2016). Vertical 
half-lines w = ± 7 r + iu, 0 < u < 00 are mapped into a branch cut i < C < ioo. In a 
similar way, vertical half-lines w = ± 7 r + iw, —00 < u < 0 are mapped into a branch cut 


—ioo < C < —i. However, 27r—periodicity of z{w) (1.5) allows to ignore these two branch 
cuts because z{w) is continuous across them. Gomplex infinities w = ±ioo are mapped 
into C = ±i. An unbounded interval [iUc,ioo), Uc > 0 is mapped into a finite interval 
[iXc,i) with 


Xc = tanh-A. 


(1.7) 


The transformation (1.6) takes care of 27r—periodicity of Stokes wave so that the function 


z(C) defined in the complex plane ^ G C corresponds to the function z{w) defined in the 
strip —TT < Re{w) = u < tt. Here and below we abuse notation and use the same symbol 
z for both functions of C, and w (and similar for other symbols). The additional advantage 
of using the mapping ( |1.6| ) is the compactness of the interval (ixc, i) as mapped from the 
infinite interval (iUc,ioo). Note that the mapping (1.6) is different from the commonly 
used (see e.g. Schwartz (1974); Tanveer (1991); Williams (1981)) mapping C, = exp (—iw) 
(maps the strip — tt ^ Re{w) < tt into the unit circle). That exponential map leaves the 
interval (iVc,ioo) infinite in f plane. 

The main result of this paper is that it was found an infinite number of sheets of 
Riemann surface with square root branch points located at C = ±ixc starting from the 
second sheet (the first sheet has the singularity only at C = iXc)- At each sheet (except 
the first one) these singularities are coupled through complex conjugated terms which 
appear in the equation for Stokes wave. In contrast, the only singularity at (j = iXc of the 
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first sheet (besides the singularity at C = i) does not have a complex conjugated sister 
at = —iXc which makes that (physical) sheet distinct from all others. It is conjectured 
that the leading order form of non-limiting Stokes wave has the form of the infinite 
number of nested square root singularities. These nested square roots form 2/3 power 
law singularity of the limiting Stokes wave as Xc 0. 

The paper is organized as follows. In Sectiona closed nonlinear integral equation for 
Stokes wave in terms of the density (jump) at the branch cut is derived and the numerical 
method to solve that integral equation is given. Section [^provides an alternative form for 
the equation of Stokes wave. Sectionj^uses that alternative form to find an asymptotic of 
both Stokes wave at Im(w) —>■ -l-oo and the jump at the branch cut. Section discusses 
a numerical procedure to analyze the structure of sheets of Riemann surface for Stokes 
wave by the integration of the corresponding nonlinear ordinary differential equation 
(ODE) in the complex plane. Section derives the analytical expressions for coupled 
series expansions at </ = ±ixc to reveal the structure of Riemann surface for Stokes 
wave. Section [T] analyzes possible singularities of Stokes in all sheets of Riemann surface 
and concludes that the only possible singularity for finite value of w is the square root 
branch point. Sectionprovides a conjecture on recovering of 2/3 power law of limiting 
Stokes wave from an infinite number of nested square root singularities of non-limiting 
Stokes wave in the limit Xc —>■ 0. In Section|^the main results of the paper are discussed. 
Appendix |A| shows the equivalence of two forms of equation for Stokes wave used in the 
main text. Appendix [B] relates different forms of equation for Stokes wave in the rest 
frame and in the moving frame. Appendix [^provides tables of the numerical parameters 
of Stokes wave. 


2 . 


Closed integral equation for Stokes wave through the density at 
the branch cut 


The equation for Stokes wave was derived in Ref. Zakharov & Dyachenkov (1996) and 
Part I (Dyachenko et al. ]M6l ) from Euler’s equations for the potential flow of ideal fluid 
with free surface (see also Appendices]^ and [^. That equation is defined at the real line 
w = u and takes the following form 


- c^Vu + yVu + H[y{l + Xu)] = 0, (2.1) 

where 

1 f(u') 

HHU) = -p.v, /LXd.' (2.2) 

is the Hilbert transform with p.v. meaning a Cauchy principal value of integral and 
subscripts in t and u mean partial derivatives here and further. The Hilbert operator H 
transforms into the multiplication operator 


iHf)k = isign(fc)/fe, 


(2.3) 


for the Fourier coefficients (harmonics) fk, 


fk = 


1 

27 T 


n 

j f{u) exp {—iku) dw, 


— 7T 


(2.4) 
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of the periodic function f{u) = f{u + 2n) represented through the Fourier series 


/(w) = ^ fk exp (iku). 


(2.5) 


k— — c 


( 2 . 6 ) 


Here sign(fc) = —1, 0,1 for /c < 0, k = 0 and fc > 0, respectively. Equation (2.5) implies 
that 

ffV = -(/-/o), 

where /o is the zeroth Fourier harmonic of /. 

It is convenient to decompose the Fourier series (2.5) as follows 

f{u) = /+(u) + f~{u) + fo, 

where 


f+{w) = ^ /fe exp (ikw) 

k^l 

is the analytical function in C“^ and 

-1 

f~{w)= y] fk exp (ikw) 


(2.7) 

( 2 . 8 ) 

(2.9) 


k— — c 


is the analytical function in the lower complex half-plane C . Then equation (2.3) implies 
that 

Hf = i{f+-r). 


Also using equation (2.3) we define the operator, 

P=^(l + ii7), 


( 2 . 10 ) 

( 2 . 11 ) 


projecting any 27T-periodic function / into a function which has analytical continuation 
from the real line w = u into C“ as follows 


pf = r + ^- 


One can apply H to (2.1) to obtain the following closed expression for y, 

[c^k - i) y - f ^ + yky^ = 0, 


( 2 . 12 ) 


(2.13) 


where k = —duH = V—and we used the following relations 

yu = Hxu and = -Hy^, (2.14) 

which are valid for the analytic function 5„(w) satisfying the decaying condition Zu{w) —)■ 
0 as u —>■ — 00 . We also assume in deriving equation (2.13) from equation (2.1) that 

(2.15) 


,Mdx = /„(„K(.)d., = o, 


meaning that the mean elevation of the free surface is set to zero. Equation (2.15) reflects 


a conservation of the total mass of fluid. Equation (2.13) was derived in Ref. Babenko 


(1987 

) and later was independently obtained from results of Ref. 

Dyachenko et al. 

(1996 
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similar equation. Ref. Babenko 


in Dyachenko et al. (2013a). See also Ref. Zakharov & Dyachenkov (1996) for somewhat 


(19871 and subsequent developments in Refs. Buffoni 


et al. (2000); Buffoni & Toland (2001); Plotnikov (1991); Shargorodsky & Toland (2008) 


used equation of the type (2.13) for the analysis of bifurcations. 


Equation (2.13) is convenient for numerical simulation of Stokes wave because it de¬ 


pends on y only as detailed in Part I (Dyachenko et al. 2016). The operator k is the 
multiplication operator in Fourier domain which is straightforward to evaluate numeri¬ 
cally using Fast Fourier Transform. 

In this paper it is however more convenient for analytical study to rewrite equation for 
Stokes wave in terms of the complex variable z. For that we apply the projector operator 


P (2.11) to equation (2.1) which results in 


c z.,, 


= -'^P [{z - S){1 + Zu)], 


(2.16) 


where f{u) = / means complex conjugation of the function f{u). Note that the complex 
conjugation f{w) of f{w) in this paper is understood as applied with the assumption 
that f{w) is the complex-valued function of the real argument w even if w takes the 
complex values so that 

f{w) = f{w). (2.17) 

That definition ensures the analytical continuation of f{w) from the real axis w = u into 
the complex plane of w S C and similar for functions of C € C. If the function f{w) is 
analytic in C“ then J{w) is analytic in C"*" as also follows from equations 

A numerical convergence of Pade approximation to the continuous density p{x) of the 
branch cut was shown in Part I (Dyachenko et al. 2016) together with the parametrization 


of the branch cut of Stokes wave as follows 


i 

5(0 = ij/6 + J 


Pix')dx' 

C-ix' ’ 


(2.18) 


where yb = yiu)\u-±Tx G K is the minimum height of Stokes wave as a function of x (or 
in the similar way as the function of u). The density p{x) is related to the jump IXjump 
of 5(C) for crossing the branch cut at C = ix in counterclockwise direction as follows 


^jump — ^(C) 5(C) — 27rp(y;), (2.19) 

see also Part I ({Dyachenko et ai 2016) for more details on that. We now use the 
parametrization (2.18) to study the Stokes wave equation (2.16). We eliminate the con¬ 


stant iyb at C = oo in (2.18) by introducing a new function 

1 


f{u) = z{u) -iyb = J 

Xc 

together with the complex conjugate 


Pix')dx' 
C-ix' ' 


( 2 . 20 ) 


/(m) = 


p{x!)dx! 

C + ix' 


( 2 . 21 ) 


which was evaluated using the definition (2.17). 


Equation (2.16) in the new valuable (2.20) takes the following form 
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- ic^fu + iyt + 2iyfc/„ + //„ + Pf-Pf-P [ffu] = 0 

with / and / given by equations (2.20) and (2.21), respectively. 

2.1. Projection in f plane 


( 2 . 22 ) 


The projector P (2.11) is defined in terms of the independent variable u. Using equation 
(2.22) together with the definition (2.20) suggests to switch from u into the independent 
variable f. To identify how to compute P in complex C-plane, we start from the Fourier 
series (2.4), (2.7) in variable u and make a change of variable (1.6) (assuming that 
—TT ^ u ^ TT and C G K) as follows 

OO OO OO 

/(u) =/(C) = = X! /nexp[2ifcarctanC] = ^ fk 


k— — oo 


k— — oo 


C-i 
C + i 


(-l)^ 

(2.23) 


where we abuse notation by assuming that /(C) = fiu) and removing ~ sign. Equations 
(2.11), (2.12) and (2.23) imply that P removes all Fourier harmonics with positive n and 
replaces the zeroth harmonic fo by /o/2 as follows 

OO p —1 

Pf{u)= = y + X! /fe exp [iA:2 arctan C] 


= V+ E A 


k——o 

-1 


n——oo 


C-i 
C + i 


(- 1 )^ 


(2.24) 


Consider a particular case f{u) = X G R and xP 0. We calculate fk by equation 


(2.4) and (2.23) through the change of variable (1.6) implying du = ^P^df as follows 


f-k = ^ f{u)e^^^du = — 


27r 


27r /_oo C - iX VC + i 


C-i 


(- 1 )* 


C" + l 


df. 


(2.25) 


Assuming k ^ 0 and closing the complex integration contour in the upper half-plane of 
C we obtain that 


f-k = i 


X- 1 
X+ 1 


(- 1 )^ 


-X + 1 


^(x) + 4,0: 


I-IX 


For the zeroth harmonic /o, equation (2.26) results in 

isign(x) 


fo — 


1 + Xsign(x)’ 


where sign(x) = 1 for x > 0 and sign(x) = —1 for x < 0- 
Using now equations (2.24), (2.26) and (2.27) we find that 


1 


-isign(x) 


C-iX 2[1-hxsign(x)] 


C + i 
C-i 


OO 

E 

k=0 
\k 


X- 1 
X + 1 


(-1)'=^-t4x)+4.o^ 
-X +1 1-1 


IX 


1 


1 


(-1)'^ = 7—^ix) + 4-x) - 


i sign(x) 


(2.26) 


(2.27) 


C - iX 


i-iX 


2[1 + xsign(x)] 


, (2.28) 


where 6{x) = 1 for x > 0 and 6{x) = 0 for x < 0. 
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In a similar way, for /(C) = 


1 


1 


(C - ix) 


2 we find from the series (2.23) that 


f-k = tt 


1 /c-i 


271-/_oo (C-ix)^ VC + i 


(- 1 )* 


c^ + i 


dC, fc > 0. 


(2.29) 


equation (2.29) that 


Closing the complex integration contour in the upper half-plane of / one obtains from 

(2.30) 


f-k = i- 


. d / C ~ i 


. \ k 


rfC VC + i 


and 




fo — - 


Jk,0 


C=iX 


1 


(i - ix)' 


r, k^O 


[1 -f Xsign(x)]^' 


(2.31) 


Taking a sum over k in equation (2.24), using equations (2.30) and (2.31) and we find 
that 


P 


1 


1 


r6'(x) 




1 


(C-ix)^ (C-ix)^“''^' ' (i-ix)^"' ' 2[l-kxsign(x)]^' 

2.2. Integral representation of the equation for Stokes wave 


(2.32) 


Using equations (2.20), (2.28) and (2.32) we obtain the following projections in terms of 

P(x): 


Pf = 


1 1 
f P{x)dx _ f ^P{x)dx 
I C-*X y 2(l-kx)’ 


Pf = - 


W{x)dx 

2(1 + X)' 


(2.33) 


We now find P [//«] used in (2.22). Equation (1.6) results in the following expression 


1 1 


c" +1 7. e + i t f p{x')p{xndx'dx'' 

j ju — 2 y *1C — 2 


(C + ix')(C-ix")^' 


(2.34) 


We perform the partial fraction decomposition of the integrand of (2.34) as follows 

C^ +1 


1 1 - x'- 


1 -1 - 2x'x" - X"" 


2(C + ixO(C-ix")" C + ix'2(x' + x")2 C - ix" 2(x'+ x")' 

1 i(i - x!'^) 
(C - ix")" 2(x' + X") ■ 


(2.35) 


and apply the projector P to (2.35) which gives with the use of (2.28) and (2.32) the 
following expression: 


P [ffu] = 


1 1 - x'- 


-k 


1 


i + ix'4(x' + x")" VC-ix" 2(1+ x"); 2(x'+ x") 


-1 - 2x'x" - X"" 


1 


(C-ix")" 2(l + x")V 2(x' + X")J 


i(l - x"") 


P{x')p{x")dx'dx" ■ 

(2.36) 
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P{x')dx' f p{x”)dx" 




(C - ix")^ ’ 


The constant yt is determined from equation (2.15) as follows 

TT OO _ 

(/-/) 


y{l + Xu)du = 


Vb 


2i 


1 


(1 + C^) 


ifc + fc) 


2dC 

J 1 + C" 


= 0 , 


which results using equation (2.20) in the following equation 

} } p{x')p{x'')dx'dx" 

Vb = 


2(x'+x"r 


P(x')dx' 

1 + X' ’ 


(2.37) 


(2.38) 


(2.39) 


Equation (2.39) allows to find yb from a given p{x)- This equation also provides a con¬ 
venient tool to estimate the accuracy of recovering p{x) by Fade approximation. For 
that one compares the numerical value of yb obtained from the Stokes solution in Part I 
( Dyachenko et al.|2016 ) with the result of the direct numerical calculation of right-hand 
side (r.h.s.) of equation (2.39) with p{x) obtained from Fade approximation in Part I (all 
these numerical values are given in tables of Part I (Dyachenko et al.|2016), through the 


electronic attachment to Ref. Dyachenko et al. (2015a|) and at the web link Dyachenko 


et al. (20156)). 


Integrating equation (2.20) in u over 27r-period one obtains the zero Fourier harmonic 
j/o of y{u) as follows 

P{x')dx' 


Vo = yb + 


1 


■X 


(2.40) 


Requiring that equations ( |2.20[ )-( |2.22[ ), ( |2.33[ ), ( |2.36[ ), ( |2.37[ ) and ( |2.39[ ) are satisfied for 
—OO < f < OO we obtain a system of equations to find the density p{x) along the branch 
cut for each c. That system has a form of nonlinear integral equation for the unknown 
function p{x)- Taking the limit C oo in that system results in the following compact 
expression 

1 r 1 1 1 


y /p(xW + 2 -J J 


Pixlpixldx'dx" 

2(x'+ x'r 


1 - 


p{x"')dx'' 


Pjxldx' 

1 + X' 


+ 


Pix')dx' 

1 + X' 


(1 - X')pix')pix")dx'dx" 

2"(x' + x")2 


= 0, (2.41) 


Xc Xc 


which can be used to find c from the given p{x)- 

2.3. Numerical solution for Stokes wave based on the integral representation 

To solve the system (2.20)-(2.22), (2.33), (2.36), (2.37) and (2.39) numerically we use the 
approximation of the integral in equation ( 2.18[ ) by the following numerical quadrature 

1 


rr \ ~/r\ ■ f Pix')dx' 

f{u) = z{Q -m = J ^ - 2^ 




C - ix' ^ c - iXi 


(2.42) 
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which has a form of Fade approximation at the discrete set of points Xc < Xi < X 2 < 
■ ■ • < Xn < 1 with weights 7 ^, j = 1,2,... N. Then the analysis of Sections 2.1 and 2.2 


with equation (2.18) replaced by the approximation (2.42) is carried out in exactly the 


same way as in equations (2.20)-(2.39) with each time p{x)dx x replaced by 7 j and 
Xj, respectively. Also integrals are replaced by the summations. It results in the discrete 
versions of these equations including 


N 




i=i 


7j 

C-*x 


N 

-E 

i=i 


2(1+ X)’ 


1 ^ 

ffu = -Ai+e)Y. 


13 ' 


N 

pf= V 
^ ^^ 2 (l + x)’ 

(2.43) 

N 


Xj" 

(C - iXT')^ 

(2.44) 


and 


N N 

p[ffu] = Y. E 

,-l-2xyXi"-X^ , 


1 


1 - xl 


1 


2(Xi' +Xi")^ 


i + ixx 4(xy + Xj"Y VC - iXi" 2(1 + Xj") 

1 1 \ i(l-X^) 1 

(C - iXi")^ 2(1 + Xi")V 2(Xj' + Xj") 


(2.45) 


Also equation (2.39) is replaced in the same discrete approximation by the following 
equation 


N N 


Vb 


= -EE 




N 


2(X/ + XT' 1 + Xj 


-E 


7t 


(2.46) 


Choosing numerical values jj, Xji j = 1,2,... ,A^ from Fade approximants of Fart 
I ( jPyachenko et aL||2016|) (these appro ximants are also available through the electronic 
attachment to Ref. Dyachenko et al. (2015a|) and at the web link iDyachenko et al. 


(20156)) we checked that equation (2.22) (together with equations (2.43)-(2.46)) is valid 
for each value of H/X with the same numerical precision as the precision (at least 10“^®) 
of the Stokes solutions of Fart I ( [Dyachenko et ar||2016 ). Values of N in Fart I range 
between tenths for moderates values of i7/A up to iV = 92 for the highest Stokes wave 
considered (given by Table 4 in Fart I ( Dyachenko et a/.|2016 )). These moderate numbers 
is in sharp contrast with the large number M of Fourier modes required for constructing 
these solutions with the same precision (M ~ 1.3 • 10® for the highest Stokes wave 
considered in Table 4 of Part I). An explanation for that dramatic difference between 
required numerical values of M and N follows from Part I. It was found in Part I that the 
error of Fourier method scales as oc exp (— 2xcAf) while the error for Fade approximation 
of Part I is oc exp {—ciXc M), ci ^ 1. It suggests that solving equations (2.22), (2.43)- 


(2.46) for numerical values of 7 ^, Xj) j = 1? 2,..., V is the attractive alternative to the 
numerical methods of Part I. 


To solve equations (2.22), (2.43)-(2.46) numerically, we aim to approximately satisfy 
equation (2.22) at the discrete set of points —00 < ( = Q < 00 , i = 1,2,...,Mi. It 
results in the nonlinear algebraic system of equations to find 7 ^-, Xji J = 1 , 2 ,... ,V. 


That system is overdetermined (see e.g. Ref. jWilkening fc Yu (2012) as the example of 
using of overdetemined systems for simulating water waves) provided we choose Mi > 2N 
but it can be solved in least square sense (by minimizing the sum of squares of the left- 
hand side (l.h.s.) of equation (2.22) taken over points C = Q, i = 1,2,..., Mi). However, 
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the difficulty in such most straightforward approach is in extreme ill-conditioning of the 
resulting algebraic system mainly because of denominators containing large powers of Xj 
clearly seen if we try to bring equation (2.22) to the common denominator. We bypass that 
difficulty here by providing the explicit procedure to find the appropriate values of Xj ; j = 
1, 2,... for each Xc (see the description of that procedure below in this Section) and 
only after that we solve equations (2.22), (2.43 1 -(2.46 1 for unknowns Xj, j = 1,2,...,iV 
at the discrete set of points —oo < ( = Q < oo, i = 1,2,..., Mi. Then the resulting 
system is the cubic polynomial in Xj, j = 1,2,...A^. That system is still moderately 
ill-conditioned but that difficulty is easily overcome by choosing Mi large enough with 
Newton’s iterations used to find numerical values of 7 ^, j = 1,2,... N, thus forming 
least-square-Newton (LSN) algorithm. E.g., for H/\ — 0.1387112446... (corresponds to 
Xc = 3.0056373876... • 10“^, see also table 0 of Appendix for details on numerical 
Stokes waves) we found that it is sufficient to use Mi = 800 and TV = 51 to achieve 
10“^® accuracy for Stokes wave. For steeper Stokes waves with H/X = 0.1401109676... 
(Xc = 6.99513864872 .. .-W^) and H/X = 0.1408682599... (xc = 5.6590609636 .. .TO-^) 
we found that using Mi = 1600, TV = 61 and Mi = lO"*, N = 78 allow to achieve 10“^® 
and 10“^® accuracy, respectively. Here values of N were chosen the same as for the 
respective Stokes wave in Part I while Mi is by a factor ~ 80 smaller than M = 65536 in 
the first case and by a factor ~ 200 smaller than M = 2097152 in the third case (values 
of M are given in Part I, through the electronic attachment to Ref. [Dyachenko et aL\ 
(2015a) and at the web link Dyachenko et al. (20156)). In these examples, using the 
symmetry of Stokes wave, the points Qj were chosen to have nonnegative values with the 
first 300 points uniformly spaced as C,i = if — 1)271/M, i = 1,..., 300 and the remaining 
Ml — 300 points uniformly (in u) spanning the remaining interval of positive values of f. 


After values of 7 ^, j = 1,2,... N are found from LSN algorithm, equation (2.421 provides 
Pade approximation for Stokes wave at the entire real line of f. Then one can use the 
results of Section 6.1 to find the high precision numerical approximation of Xc which 
completes the current step in 77/A (or equivalently the current step in Xc)- These step 
are repeated to gradually increase H/X (or equivalently decrease Xc) by changing the 
velocity parameter c to span the desired range of Stokes waves. 


The procedure to find the grid Xjj j = 1,2,... N a.t each step is the following. Assume 
that Xc < Xi < X 2 < • • -Xn-i < Xn < 1 and Xc ^ 1- We use the property of Stokes 
wave that p(x) changes a little vs. a change of i7/A for x Xc provided Xc ^ 1- It 
implies that Xj can be chosen independently on Xc for all j such that Xj ^ Xc- In 
numerical examples above we chose numerical values in the range Xj ^ Xc from Pade 


data for Stokes wave with H/X = 0.1409700957... obtained in Part I (Dyachenko et al. 


2016). Also the grid Xj can be chosen from the grid obtained at the previous step (with 
a previous smaller value of Xc)- 

We now consider the construction of grid for smaller values of Xj- H assume a 
power law singularity p{x) oc (x~ Xc)“j a > 0 and consider the limit (/ — >■ ixc in equation 
(2.18), then the transformation to a new integration variable t = (x — Xc)"^ 


removes 


the singularity from the integrand in equation (2.18). The uniform grid tj = jAt, j = 
1,2,..., At = const in t is the natural choice to use for the integration in the variable t. 
The corresponding grid in x is given by 


X,-Xc = iy“=//“AtV^ 


(2.47) 
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Stokes wave has the square root singularity at ^ = ixc with the expansion 

OO 

5(C) - m = /(C) = E - iXc)^■/^ 


(2.48) 


j=o 


where aj are real constants (see Part I ( Dyachenko et aL|2016 ) as well as Sectionsand 
below for the justification of that expansion). It implies (see Part I) the square root 
singularity for the density p{x) oc {x~ Xc)^'^^ in the integrand of equation ( 2.18[ ). Using 
equation (2.47) with a = 1/2 one then obtains that 


Xj -Xc = j -- 1, 


(2.49) 


which is in the excellent agreement with numerical values of Xj obtained in Part I provided 
^ O.OlXc- 

In the range Xc ^ X 1j the density p{x) oc x^^^ is well approximated by the density 
p{x) ot x^^^ of the limiting Stokes wave (1.1) as shown in Figure 8 of Part I. Using 
equation (2.47) with a = 2/3, one obtains that 


X, = (2.50) 

where c is the positive constant and j ^ 1 such that Xj 1 - We additionally have to 
approximate the transition between two scalings (2.49) and ( |2.50 ) at the intermediate 
values of j. Exploring fits of Xj '^s. j for multiple sets of numerical data of Part I we 
found that a satisfactory fit (including the required transition) is given by the linear 
combination of the scaling (2.49) and ( |2.50 1 superimposed with the exponential growth 
in j as follows 

Xj = Xc Cif + , (2.51) 

where the positive fitting constant Ci, C 2 and C 3 changes slowly with Xc (change of Xc in 
5 orders of magnitude results in change of these constants by less than 50%). 

Based on these observations we implemented the following procedure to find numerical 
values of ci, C 2 and C 3 for each value of Xc- We choose N from the previous step (with 
the previous value of Xc)- (If performing the current step we are not able to reach the 
desired precision with the increase of Mi, i.e. LSN algorithm would not converge to the 
prescribed tolerance, e.g. 10“^®, then N has to be increased by 1 and the current step 
restarted from the beginning). Next we choose jmatch from the values Xj of Table 4 of 
Part I (or from the grid obtained at the previous step in Xc) such that Xj^atahlXc ^ 100 
which well ensures the required condition Xj ^ Xc- (For larger values of Xc one can 
use a smaller value of Xjmo.tah.IXc to make sure that Xjmotah 1 - F.g. for the case 
H/X = 0.1387112446... (first numerical case mentioned above in this section) we choose 
XjmotahIXc ~ 34.) After choosing the value of j match j we perform 4th order interpolation 
of Xj Ets the function of j and find values of the first and second derivatives, x'j s-nd 
Xj, of that interpolant at j = jmatch- We use these 3 numerical values Xjmotch^ 
and x'j.^„f^h fo find the numerical values of ci, C 2 and C 3 by matching the corresponding 


values of equation (2.51) and its two derivatives at j = jmatch- Then equation (2.51) 


provides the numerical values of Xji J = 1)2, ... , jmatch completing the construction of 
the numerical grid Xj, / = 1 ) 2 ,. ■., for the current step in Xc- (If any of the constants 
ci,C 2 or C 3 turns negative then one has to decrease XjmotahfXc to avoid that but in our 
numerical examples we experienced such problems only if Xjmotah/Xc was chosen > 10 ^). 
Then LSN algorithm is used as described above. 

The efficiency of the grid Xj, j = 1,2,... A thus obtained requires a good initial 
estimate of Xc with the relative accuracy ~ 10“^. It is achieved by a gradual increase 
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of H/X (decrease of Xc) at multiple previous steps of LSN algorithm. Values of Xc are 
found at each previous step with high precision by the procedure of Section |6.1[ The 
polynomial extrapolation of Xc to the current step is performed to reach the needed 
relative accuracy ~ 10“^. Note that the numerical detection of the incorrect value of Xc 
prediction is straightforward because it would cause the oscillations of (with changing 
of its sign) around several smallest values of the index j = 1,2,3,.... 

We would like to stress the difference of LSN algorithm of this Section compare with the 
method of Part I ( Dyachenko et aL|2016 ). Stokes wave was obtained in Part I by using the 
Fourier series representation of the solution combined with Newton-Conjugate-Gradient 
iterations method. After that the resulting solution was approximated at the real line 
of ^ by a numerically stable version of the Pade algorithm. Thus Pade approximants 
of Part I were only the auxiliary tool to compactly represent the result of calculation 
of Stokes waves. In contrast, in this Section we completely bypass the Fourier series 


representation and numerically solve integral equations (2.22), (2.43)-(2.46) directly by 


Pade approximants. The cost of the approach of this Section is that instead of 2N free 
parameters Xjjlj^ V = 1, 2,..., V of Pade approximats of Part I, we now have only N 
free parameters y^-, V = 1, 2 ,..., V, while the values ofxj, V = 1, 2 ,..., V are fixed by 
the grid algorithm described previously in this section. It means that to achieve the same 
precision we need to approximately double the value of N compare with Part I. This is 
however very moderate cost compare with the Fourier method of Part I. 

In conclusion, in this Section we demonstrated the performance LSN algorithm for 
several values of Xc which were previously explored in Part I by the Fourier method. We 
expect that much smaller values of N, required for LSN algorithm compare with Fourier 
method, will allow to find Stokes waves for much smaller values of Xc than achievable by 


Fourier method of Part 1. In addition, equation (2.41) can be used to exclude c from the 


system allowing to gradually increase H/\ (or equivalently decrease Xc) thus avoiding the 
problem of nonmonotonic dependence of c on H/X encountered in Part 1. The detailed 
practical realization of that limit of smaller Xc is beyond the scope of this paper. 


3. Alternative form for the equation of Stokes wave 

The equation for Stokes wave can be written in a form which is alternative to equation 


(2.1) as follows 


2/ = - 2) = Y ( 1 - 


1 


(3.1) 


Appendix [A| shows the equivalence of both forms of equations (2.11 and (3.1) for Stokes 
wave. Also Appendix [b] discusses differences in derivation of equations (2.1) and (3.1) 
from basic equations of the potential flow of ideal fluid with free surface. Different versions 


of equation (3.1) (up to trivial scaling of parameters and shift of y by different constants) 
were used by Grant (1973), Williams (1981), Plotnikov (1982) and Tanveer] ( 1991| ). 
Transforming equation (3.1) into the variable f (1.6) results in 


Z — Z = 1C 


1 - 


Solving equation (3.1) for z, 


(l+C^)2|zcP 
one obtains that 


(3.2) 


z„[i(z - z) + c2] 


(3.3) 


which is the nonlinear ODE provided z„ is known. In a similar way, equation (3.2) results 
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= 


(1 + Zc^[\{z - z)+ c^]' 


(3.4) 


Equations (3.3) and (3.4) can be considered as ODEs for z(u) and z(C), respectively, 
if z is the known function. Then solving ODE provides a convenient tool to study the 
analytical properties of Stokes wave in different sheets of Riemann surface of z. 


4. Asymptotic of Stokes wave at Im{w) —?■ +cx) and jump at branch cut 

An asymptotical solution of Stokes wave in the limit Im{w) —> +oo is obtained from 


equation (3.3) as follows. Equation (2.9) implies the exponential convergence oc e of 


z{w) to its zeroth Fourier harmonic, z{w) ij/o for Im{w) —>■ —oo. Here yo is determined 
by the mean-zero elevation condition ( |2.15 ) and is given by equations (2.39) and (2.40). 
Respectively, z(w) converges exponentially to —ij/o for Im{w) —?> oo. Then Zu and z in 
equation (|3.3|) can be replaced by 1 and —ij/Oj respectively in that limit resulting in 


1 + — 


i(z -b iyo) + 


Im(w) 


oo. 


(4.1) 


Integrating equation (4.1 1 in the upper right quadrant w G C"*", Re{w) > 0 for u ^ 1, 
one obtains that 

z{w) — ic^ In [z(w) -I- iyo] = —w + c+, (4.2) 

where c_|_ is the constant. A similar integration in the upper left quadrant w G C"*", Re(w) < 
0 , for u >> 1 results in 

z{w) — i(? In [z{w) -|- ij/o] = —w + c_, (4.3) 

where c_ is the constant. 


Taking w = tt -|- iu in equation (4.2) and w = —tt -|- iu in equation (4.3) together with 


the periodicity condition z(7r -|- iu) = z(— tt -|- iu) result in the condition for constants c+ 
and c_ as follows 

c_ — c+ = —27T. 


(4.4) 

Exponents of equations (4.2) and (4.3) are similar to the Lambert IT-function. Solving 


these equations in the limit u —>■ oo (see e.g. Refs. Dyachenko et al. (20136); Lushnikov 


et al. (2013) for details on a similar technique) one obtains that 


z{w) = 


, , • 2, r , , • 1 + c± -b iyo] 

— w -b C± -b 1C In [—w -b c± -b lyo-^- 

-w + c±+ lyo 

ln[-u;-bc± -biyo]'' ^ 


+0 


-w -b C± -b iyo 


(4.5) 


where a use of c+ and c_ assumes that Re{w) > 0 and Re{w) < 0, respectively. If 
equation (3.3) is used instead of the reduced equation ( |4.1[ ) in derivation of equation 
(4.5), then an additional exponentially small error term 0{e~^/v) appears in r.h.s. of 
equation (4.5). The two leading order terms —w + c± and ic^ In [—rc -b c± -b iyo] in r.h.s. 


of equation (4.5) are similar to equation (2.22) of Ref. Tanveer (1991), where these terms 
were derived in somewhat similar procedure to the derivation of equation (|4.5|). 


One concludes from equation (4.5) that z{w) has a complex singularity at z = oo which 


involves logarithms with the infinite number of sheets of Riemann surface. Full analysis 


of that singularity requires to study next order terms in equation (4.5) which is beyond 
the scope of this paper. 
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Taking an additional limit w = iv ± e, e > 0, e —?► 0 in equation (4.1), using the 


condition (4.4) and expanding in u ^ 1, one obtains the jump at the branch cut 


z(w — 0) — z{iv + 0) = —2n + 


2nc^ 


+ 0 (u-") 


(4.6) 


where the branch cut v S [ivc,ioo] is crossed in counterclockwise direction. 

According to equation (2.19), the jump (4.6) is related to the density p{x) ( |2.18 ) as 
follows 

c2 

P(x)lx=tanh(«/2) = 1-^), ■U>1. (4.7) 


For u » 1, one obtains from equation (1.6) that 1 — y 1 and u = — In +0(1 — 
x). Then the density (4.7) takes the following form 


p{x) = 1 + 


O 


Vln''(l-x) 


1 


Equation (4.8) implies a unit value 


p(l) = 1 


and the divergence of the derivative 

dp{x) 


dx 


(I-X)ln^(ii^) 


oo for X —t 1- 


(4.8) 

(4.9) 

(4.10) 


5. Numerical procedure to analyze the structure of sheets of 
Riemann surface for Stokes wave by ODE integration 

We use Fade approximants of Stokes wave found in Part I ( jPyachenko et aL|2016 ) and 
provided both in tables of Part I, through the electronic attachment to Ref. DyachenT^ 
et al. (2015a I and at the web link Dyachenko et al. (20156) in the following form 


N 


2^parfe(^) — U + iyi} + ^ ( 


^ tan(u/2) - ixr, 


(5.1) 


with the numerical values of yb, the pole positions Xn and the complex residues 7 „ 
(n = 1,..., N) given there. These data of Pade approximation allow to recover the 
Stokes wave at the real axis w = u (and similar at C = Re{f) in the complex C-plane) 
with the relative accuracy of at least 10“^® (for the vast majority of numerical cases the 
actual accuracy is even higher by several orders of magnitude). 

Analytical continuation of the Pade approximant ( |5.1| ) from utowGCis given by the 
straightforward replacing of u by w. That analytical continuation is accurate for w S C“ 
but looses precision for w € C’*' in the neighbourhood of the branch cut w € [iUc,ioo) 


where the discrete sum (5.1) fails to approximate the continuous paramterization (2.18) 
of the branch cut. Thus a significant loss of precision compare to 10“^® occurs only if 
the distance from the given value of f to the branch cut is smaller or comparable with 
the distance between neighbouring values of Xn in equation (5.1). 


Numerical integrations of ODE (3.4) (and occasionally ODE (3.3)) in this Section were 


performed using 9(8)th order explicit Runge-Kutta algorithm with adaptive stepping 
embedded into Mathematica 10.2 software. That algorithm is the implementation of 


Ref. Verner (2010) and is based on the embedded pair of 9th and 8th order methods 


with higher order method used for the adaptive step-size control. We used the numerical 
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precision of 55 digits and reached the accuracy 10“^° to make sure that no significant 
accumulation of ODE integration error occurs in comparison with 10“^® precision of the 
Fade approximants of Stokes waves. We also independently verified the accuracy of the 
numerical ODE integration by comparing with the analytical results of Section in the 
neighborhoods of ^ = Hxc in multiple sheets of Riemann surface. 

At the first step of our investigation, ODE (3.4) was solved numerically to find the 
approximation Zode{C) for z{Q with C G C’*' in the first and the second sheet of Riemann 
surface using the approximants (5.1) for z and Z(^. Here the first (physical) sheet of 
Riemann surface corresponds to z{C) with fluid occupying C G C“. The second (non¬ 
physical) sheet is reached when the branch cut C G [iXcj i] (or equivalently w G [iuc, ioo)) 
is crossed from the first sheet. That ODE was solved with initial conditions at real 
line C = i?e(C) by integrating along different contours in C S C"*". A high precision of 
at least 10“®° was achieved in ODE solver to avoid any significant additional loss of 
precision compare with 10“^® precision of equation (5.1 1 . That ODE solution used Zpade 
and {zpade)c which through the complex conjugation corresponds to the approximants 
(5.1) in C € thus avoiding any loss of precision compare with 10“^®. We stress 
here that the use of Fade approximation is the auxiliary tool which does not make any 
difference in the final result because it matches the precision of Fourier series. The Fourier 
series of Fart I can be used directly instead of Fade approximats which however would 
require significant increase of computational resources to reach the same precision. 

The left panel of Figure]^ shows a typical OBCDAO rectangular contour for ODE in¬ 
tegration which was used for the analytical continuation of Stokes wave into the second 
sheet of Riemann surface. The ODE solution in the second sheet is obtained when in¬ 
tegrating contour crosses the branch cut [iXcO]- The second subsequent crossing of that 
branch cut returns zode{Q to the first sheet confirming the square root branch point at 
C, = ixc- Figure [^provides a numerical example of such double crossing. In other words, it 
was found that ODE integration along any closed contour in C G C"® with double crossing 
of the branch cut (twice integrating along OBCDAO ) always returns the solution to the 
original one. If the height of OBCDAO contour is made smaller than Xc then there is no 
crossing of the branch cut and OBCDAO integration returns to the initial value after a 
single round trip as shown by dashed curves of Figure In similar way, if the height 
of OBCDAO exceeds 1 then there is no crossing and zode{C) stays in the first sheet 
(crossing of the branch cut (i, ioo) corresponds to the jump on 2tt in u direction in w 
plane while there is no jump in z because of 27r-periodicity). 

We also verified that there are no singularities in the limit |i?e(C)| —>■ oo by switching 
to ODE integration (3.3) in w variable. In that limit Re{w) —>■ ±7r which allows to extend 
the contour in w over the entire 27r period in u direction (in Q variable it would require 
to integrate over the infinite interval —oo < Re{C,) < oo). For all subsequent cases in 
this section it is assumed that such integration in w was performed to check the limit 
Re{w) ±7r. 

The second step of our investigation was to find z(C) by integrating ODE (3.4) in the 
second sheet with C, G C“ using the complex conjugate of zode{0^ C £ 'C''' found at 
previous step to approximate z and . Initial condition at that step was at the real line 
C = Re(C) with z(C) obtained at the step one for the second sheet. 

The second step reveals a new square root singularity at C = —iXc in the second sheet. 
Similar to the step one, the double integration over the contour ABEFA shows that 
z(C) returns to its original value confirming that C = —iXc is the square root branch 
point. Crossing of the branch cut [—iXc ~i] (corresponds to that new branch point ^ = 
—iXc) allows to go into the third sheet of Riemann surface. At that crossing one has to 
simultaneously cross from the first to the second sheets for z and z^ which again are the 
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Figure 2. A sche mati c of integrating contours in different sheets of Riemann surface in the 
complex variable ( (1.6 1 near the origin. The fir st (p hysical) sheet has a square root singularity 
only at = ixc in (UU Then integrating ODE (3.41 over the closed contour ABC DA provides 
the analytical continuation into the second sheet of Riemann surface as the branch cut (dashed 
line) is crossed. A s a r esult, 2 (<C) does not return to its initial value at the origin 0. In contrast, 
integrating ODE (3.41 over the closed contour ABEFA (or over ABC DA provided its height 
falls below f = i^c), one does not cross the branch cut so z{C,) returns to the same value at 
the origin with 2 (^) remaining in the first sheet. The second sheet has the second square root 
branch point singularity at ^ = —iXc at the lower complex half-plane C”. Integrating over 
contour ABEFA in the second sheet results in the analytical continuation of z{Q into the 
third sheet of Riemann surface. Starting from the third sheet, extra square root branch points 
appear away from the imaginary axis. Branch cuts for these off-axis singularities are chosen to 
be extended horizontally as shown by dashed lines in the two right panels. The number of these 
branch points grows with the growth of the sheet number as schematically shown in the right 
panel. We avoid crossing these branch by modifying contours ABCDA and ABEFA as shown 
in the two right panels. Note that these two contours must by symmetric with respect to the 
real line even if the chosen pair of off-axis singularities (symmetric with respect to the imaginary 
axis) are located only in one of the complex half-p lane s C”*" and C“ in the given sheet. This is 
because 2(0 is needed for the integration of ODE (3.4|. 


complex conjugate of zode{C), C ^ found at previous step. It was found that the 
third sheet has branch points both at (^ = iXc and f = —iXc- In a similar way to previous 
steps, at the step three one crosses the branch cut [ixc,i] to go into the fourth sheet of 
Riemann surface which found to has branch points both at C = iXc and f = —iXc- At 
the step four one crosses the branch cut [—iXc —i] to go into the firth sheet of Riemann 
surface which again has branch points both at C = ixc and C, = —iXc etc. At each 
sheet, used in integration of ODE (3.4) is behind by one in sheet number to the 
current sheet, i.e. values of 2, 2^ from the first, second, third etc. sheets are used for ODE 
integration in the second, third, fourth etc. sheets, respectively. After exploring several 
hundreds of sheets for different values of Xc one concludes that the number of sheets 
is infinite. The double integration over the contours ABCDA and ABEFA shows that 
f = ixc and C = ~iXc are square root branch points in all sheets. In the next section this 
conjecture is strengthened by the analysis of expansions ai f = ±ixc in these multiple 
sheets. 

Starting from the third sheet, extra square root branch points appear away from the 
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Figure 3. The amplitude |5(C)| (a) and the argument Arg(5(^)) (b) vs. arclength in the 
variable ( (scaled by Xc) along the closed contour OBCDAO shown on the left panel of Figure 

t (the contour is passed twice in the counterclockwise direction) for ODE integration (provides 
e analytical continuation of Stokes wave in the complex plane) of Stokes wave solution with 
H/X = 0.1387112446 ... (corresponds to Xc = 3.0056373876 ... • 10~^, see Tablef^of Appendix 
Ofor details on numerical Stokes wave). The contour width is \AB\ = 2 |0B| = 2Xc. Solid lines 
are for the contour height |i3C| = 4xc (the contour OBCDAO twice crosses the branch cut 
[iXc,i] in counterclockwise direction) and dotted lines are for the contour height \BC\ = Xc/2 
(for that height the contour OBCDAO does not crosses the branch cut [ixc i] as well as the 
total arclength is smaller). It is seen that solid lines are periodic over the total arclength (two 
round trips around the contour OBCDAO are needed to return to the initial value 2(0) in the 
first sheet of Riemann surface) compare with the half-arclength periodicity of dotted lines (the 
contour OBCDAO is located in the first sheet only with one roundtrip sufficient to return to the 
initiaf vafue 2 ( 0 )). 


imaginary axis. The existence of these off-axis singularities are closely related to the 
analysis of Section]^ They are located significantly more far away from the origin than 
the on-axis = ±ixc singularities. These singularities appear at each sheet starting 
from the third one in pairs located symmetrically with respect to the imaginary axis 
as schematically shown in Figure The symmetric location of pairs of singularities are 
required from the symmetry condition 


2(-C) = -2:(C)- (5.2) 

That symmetry condition results from the symmetry y{x) = y{—x) of Stokes wave in 
physical variables. The location of the first pair of off-axis square root singularities at 
C, = Cd and C, = — Cci is schematically shown in the third panel of Figure]^ By adaptively 
increasing the horizontal and vertical sizes of the contour ABEFA of Figurej^ we found 
that for Xc “C 1 the first pair of off-axis square root singularities are located in the third 
sheet at C = Cci and C = ~Cci with 

Cci (17.1719 - i 10.7734)xc. (5.3) 


Other off-axis pairs are located even more far away both from the real and imaginary 
axes as schematically shown in Figure Branch cuts for all off-axis singularities are 
chosen to be extended horizontally as shown by dashed lines in the two right panels of 
Figurej^ For |A0|, |0i?| < 17.1719xc, one can use the same contour as in the left panel 
for all sheets. However, for larger values of |H0|, |0i3| one has to bypass off-diagonal 
singularities as shown in the two right panels of Figure to keep the enumeration of 
the sheets as described above (based on on-axis C = ±iXc singularities). The number of 
off-axis branch points grows with the increase of the sheet number. We also performed 
double integration over closed contours around multiple off-axis singularities and found 
that each of them is the square root branch point. 

By-product of ODE integration of this section is that one can also calculate the jump 
—27rp(x) (see equation (2.18)) at the branch cut of the first sheet with the high precision. 
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E.g. one can start ODE integration at C = 0 in the first sheet and integrate until reaching 
a small neighborhood of C = iXc without crossing the branch cut f G [iXoi]- After that 
one can integrate ODE independently along two line segments ( = [±e + ixc,±e + i], 
e —>■ 0 and calculate a difference between these two integrations recovering p(x) with 
the precision of our simulations ~ 10“^®. A comparison of that high precision p(x) with 
the numerical approximation of p(x) obtained in Part I from the continuous limit of 
Fade approximation (see Eigs. 6b, 7 and 8 in Part I) confirmed the numerical error order 
estimates of Section 4.2 of Part 1. Also we found that equations (4.6)-(4.10) are also in 


the excellent agreement with the numerical values of p(x) confirming the asymptotical 
analysis of Section]^ 

ODE of the type (|3.4|) was numerically integrated in Ref. Tanveer (19911 based on 


the Taylor series representation of Stokes wave in the physical sheet (the additional 
conformal mapping from the unit disk used in Ref. Tanveer (1991) into the half-plane 
C“ of C (1-6) makes that Taylor series similar to the Eourier series representation of 


Part I (Dyachenko et al. 2016)). That representation allowed Ref. Tanveer (1991) for 
the first time to extend the numerical integration into the upper half C“*" of the second 
Riemann sheet and demonstrate the existence of the square root branch point there. 


Thus a numerical ODE integration of Ref. Tanveer (1991) is similar to our first step of 
this Section restricted to C+ only. 

Note that it was assumed throughout this Section that any crossing by the ODE 
integration contour of both [i,ioo] and [—ioo,—i] is avoided. Such crossing would be 
harmless in the first sheet because of 27T-periodicity of z(w). However, starting from 
the second sheet, z{w) is generally non-periodic in w. Thus the branch cuts [i,ioo] and 
[—ioo, —i] cannot be ignored any more contrary to the case of the first sheet case discussed 
in the Introduction. It implies that a crossing of these branch cut provides the additional 
sheets of Riemann surface. We however do not explore these sheets here because they 
have the distance 1 from the real axis in f plane for any value of Xc thus not contributing 
to the formation of the limiting Stokes wave. 


6. Series expansions at (( = ±ixc and structure of Riemann surface for 
Stokes wave 


Equation (3.1) together with the definition (2.17) shows that singularities at C = ±Ix, 


are coupled through complex conjugation. We found in Part I (Dyachenko et al. 2016) 


that there is only one singularity (a square root branch point) in the first (physical) sheet 
of Riemann surface which corresponds to the finite complex w plane. In addition, there 
is a singularity at C = i which is the complex infinity w = ioo and is discussed in Section 
1^ Following Part I we chose the line segment [ixc, i] as the branch cut connecting these 
two singularities in the first sheet of Riemann surface as sketched on the left panel of 
Figure Singularity at ^ = —iXc is not allowed in the first sheet because z is analytic 
in the fluid domain w G C~. 

Consider the expansions in /th sheet of Riemann surface 


zp+iC) =^ie^■^’'/^a+,^J(C-iXc)^^^ ; = 1,2,..., 
j=o 


and 


( 6 . 1 ) 


ziAC) = E + iw)^■/^ Z = 1, 2,..., 

1=0 


( 6 . 2 ) 
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where subscripts “ + ” and “ — ” mean expansions at C = iXc and ^ = — ixc respectively. 
Here the branch cuts of (C — iXc)^^^ and (C + iXc)^^^ are assumed to extend from ( = iXc 
upwards and from ( = —iXc downwards, respectively as shown in Figure Often a 
location of the branch cut of square root is taken on the negative real axis of the argument. 
To use that standard agreement about a location of the branch cut, one can replace 
{( — and + iXcY^'^ in equations ( 6.11 and ( 6.21 by (—i)'^/^(iC + Xc)'’^^ and 

F^(—iC + XcY^'^j respectively. 

Following Section we enumerate sheets of Riemann surface according to the branch 
points ^ = ±iXc as follows. A crossing of the branch cut [iXoi] in the counterclockwise 
direction means going from I = 2n— 1th sheet of Riemann surface to I = 2nth sheet with 
n = 1, 2,.... Case I = 1 corresponds to the physical sheet of Riemann surface. Similarly, 
crossing of a branch cut [—iXc ~i] in the counterclockwise direction means going from 
I = 2nth sheet of Riemann surface to I = 2n + 1 sheet with n = 1,2,.... Plugging 
expansions (6.1) and (6.2) into equation (3.1) and collecting terms of the same order of 
(C ± iXcY^'^ result in the following relations 


n-,2n,i — 0, 

-2 

,2n,2 — Z 9 5 

l-Xc 

16c^ 

a-,2n,3 = —- ^“2 — Ti 

3 (1 - Xi) a+,2n-i,i(c^ - a+,2n-l,0 - a_,2n,o) 

2xc 4c^ 

(1-Xc^)^ (1 - Xc^)^(c^ - a+,2n-l,0 - a_,2n.o)^ 

_ 8c^[2 + (—1 + X^)a^_^2n-1,2] _ 

(1 - Xc^)^a+,2n-l,l(c^ - a+,2n-1.0 “ a-.2n.o) ’ 


(6.3) 


for n ^ 1 and 

n+,2n+l,l = 

n+,2n+l,2 = 

+ 


16 c^ 


3(1- Xc^) a_,2n.3(c2-a-,2n,0 + a+. 2 n+l,o) ’ 
2 128c'‘ 


1 - Xc^ 9(1 - Xc^)^a?. 2„,3 (c2 - a_, 2 n.O - a+, 2 n+l.o)^ 
32c^[-2Xc + (1 - Xc)^a-.2n.4] 


(6.4) 


9(1 - Xc )'‘a?. 2„_3 (c2 - a_,2n,0 - a+,2n+l.o) ’ 
n+,2n+l,3 = ■ • ■ I 


for n ^ 1. 

One cannot take n = 0 in equation (6.4) which corresponds to ^ = 1 (the physical sheet 
of Riemann surface). This special case has to be considered separately because in the 
physical sheet there is no singularity at C = —iXc (no singularity inside fluid domain). It 
implies that 


a_,p 2 i+i = 0 for j = 0,1,2,... 


(6.5) 


Solving equations (6.1), (6.2) and (3.1) for Z = 1 with the series expansion at ( = —ix, 
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2 

a+,1.0 — c — a-,1.0) 

_ _ -23/^c _ 

" (1 - [(2 + (1 - Xl)a-,i,2)f^ ’ 

_ [2 + (1 - Xc)®-.i,2)] ^ (6.6) 

2V218c(1-Xc2)3/2 

2^/^c [2xc - 2xc(-l + Xc)«-,i,2 + (-1 + xD'^a-Xi] 

(1_^^2)3/2[2+(1_^2)„_i_2)]3/2 

0+,1.4 = • • ■ , 


Expressions (6.6) are uniquely determined by values of c, Xc and a-^i^ 2 j, j = 0,1,2,..., 
where all expressions under square roots are positive and the principle branch of all square 


roots is assumed. In contrast, the expressions (6.3) and (6.4) are not the unique solutions 


of equations (6.1), (6.2) and (3.1). In addition to the solution (6.3), one can obtain two 


more spurious solutions for a__2nj - However, these spurious solutions do not correspond 
to Stokes wave. One spurious solution has a__2n,2j+i = 0 for j = 1,2,..., i.e. it does 
not have a singularity at ^ = —iXc- The second spurious solution has either a radius 
of convergence well below Xc or even the zero radius of convergence. Both solutions are 
spurious because they cannot have the same value in the region of overlap of the disks of 


convergence of both expansions (6.1) and (6.2). After spurious solutions for a-^ 2 n,j are 
discarded, one obtains the unique solution (6.3) as well as the solution (6.4) a+^2n+ij' 
also turns to be uniquely defined. Another peculiar property of the solution (6.3) is that 


a-,2n,i = 0 while a+,2ri+i,i 7^ 0 as given by the solution (6.4). 


R.h.s. of equation (6.4) provides the explicit expressions for the coefficients a_|__2n+i,jj 
j = 1 , 2 ,,..., for 2n + 1th sheet of Riemann surface through the coefficients a_^2nju 
0 ^ Ji ^ j + 2 at 2nth sheet. The only coefficient which remains unknown is the zeroth 


coefficient a+^2n+i,o for each n ^ 1. In a similar way, r.h.s. of equation (6.3) provides the 
explicit expressions for the coefficients a-,2nj j = 1,2,,..., for 2nth sheet of Riemann 
surface through the coefficients a+,2n-iji, 0 ^ ^ j — 2, j = I, 2,,..., at 2n — 1th 

sheet. The only coefficient which remains unknown is the zeroth coefficient a_,2n,o for 
each n ^ 1. 


The explicit expressions for a+_2n+i,j 


and Q-.2ri,.7 turn cumbersome with the increase 


of j beyond values shown explicitly in equations (6.3) and (6.4). The explicit expression 


a+,2n+i,j and a-^ 2 n,j were obtained with the help of symbolic computations in Mathe- 
matica 10.2 software. These expressions were used to calculate values of all coefficients 
a+,2ra+i,j and a-^ 2 n,j for J ^ 1 numerically with any desired precision (typically we used 
quadruple (quad) precision with 32 digits accuracy and took into account all j in the 
range 1 ^ j ^ 200). The remaining coefficients a_,2n,o and a+^2n+i,o for each n ^ 1 as 
well as the numerical value of Xc were determined by a numerical procedure which is 
described below in Sections 16.11 and 16.21 


Values of a+^ 2 n+ 2 ,j and a_^2n+i,j are obtained from a+_2n+i,j, and a__2n,j by the 
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(6.7) 


following relations 

^+,2n+2,j — ( 1 )^ a-|_^2n+l,j; ^ — 0; 1? ^7 • ■ • ? 

,2n+lj' ( ,2n,j ; ^ ; 

which immediately follows from the condition at the crossing of branch cuts. 

6.1. Finding of Xc^ from matching the series expansions at f = ±iXc first sheet 

Equations (6.6 1 determine values of a+q j, j = 0,1,2,. .. from a_q_ 2 j, j = 0,1, 2, ... thus 
relating the series expansions at C = —iXc and ( = ixc at the first sheet. The series at 
f = —iXc is given by equation (6.2) with I = 1 together with the condition (6.5). That 
series contains only integer powers of C + iXc- The disk of convergence |C + iyd < r of 
that series is determined by the branch point at ^ = ixc which implies that the radius of 
convergence is r = 2xc- The series at C = iXc at the first sheet is given by equations (6.1), 


( 6 . 6 ) and contains both integer and half-integer powers of C~iXc- The disk of convergence 
1C ~ iXcl < »" is determined by the branch point at C = —iXc of the second sheet. Thus 
the radius of convergence is also r = 2xc- In other words, the radius of convergence of 
the series ( 6 . 1 ), (| 6 . 6 |) in the physical sheet is determined by the singularity in the second 


(non-physical sheet). 

Numerical values of the coefficients a_,i_ 2 j) j = 0,1, 2,... are immediately obtained by 
the differentiation of the Fade approximants of Part I for each numerical value of H/X. 
Accuracy of that approximation of the coefficients a_q^ 2 j is checked by plugging these 
numerical values into the series (6.2) with I = 1 and using (6.5). For numerical evaluation 


that series is truncated into a finite sum 


,s'um(C) — 


E 

1=0 


ie-i^'"/V,i,,(C + iXc)'/" = 


E 

1=0 


ie-i^'"/^a_,i,2,(C + iXc)d (6.8) 


where jmax is chosen sufficiently large to match the numerical precision of Fade approx¬ 
imants. It is convenient to evaluate that sum at ^ = 0 which is well inside the disk of 
convergence |C + iXcl < ‘^Xc- It was found that jmax = 200 at C = 0 is well sufficient 
to reach a numerical precision about quad precision ~ 10“^^ of simulations of Part I. 
That numerical value of jmax (sufficient to reach quad precision) is only weakly depen¬ 
dent on H/X. To understand that weak dependence one can note that |C + iXclc=o = Xc 
which is one-half of the radius of convergence of the series (6.2). The asymptotics of the 


terms of the series ( 6 . 2 ) for large j is determined by the radius of convergence as follows 
\o--.i^ 2 j/o,-x^ 2 j+ 2 \ — 2xc- Then the truncation of the series (6.2) by the finite sum ( 6 . 8 ) 
= 200 gives the error ^ ~ ~ 10 “^° 


m comparison 


with jmax = ^tlU gives tne error ^ 
with Fade approximation of Stokes wave at ^ = 0. 

It worth to note here that the number of derivatives jmax/‘2 = 100 which was reliably 
recovered above from Fade approximation is really large which demonstrates the highly 
superior efficiency of Fade approximation compare with Fourier series. E.g., if instead 
Fade approximation of Stokes wave, one uses the Fourier series representation of Stokes 
wave, then the number of derivatives calculated from that series with a high numerical 
precision would be limited to just a few (about 10-20 derivatives if the relative error ~ 1 
in derivatives is allowed). 


To obtain numerical values of a+,ij-, j = 0,1, 2,... from equations ( 6 . 6 ) one also has 
to know the numerical value of Xc- Part I described a numerical procedure to recover Xc 
with the accuracy ^ 10 “^° which is significantly below the accuracy < 10 “^® of numerical 
Stokes solution itself and its Fade approximation. In this paper to greatly improve that 
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precision of Xc, one sets a condition that Xc is chosen in such a way to allow the series 
(6.1) to recover the value of z(0) with the accuracy better than ^ 10“^®. 


Similar to equation (6.8), the series (6.1) is truncated to a finite sum 


Jmax 

E 

3=0 


2i,+.s«m(C) = E - iXcF/^ 


(6.9) 


where we again choose that jmax = 200 which is well enough to match quad precision 
^ 10“^^. Contrary to equation (6.8), the sum (6.9) includes also half-integer powers 
of C ~ iXc because f = ixc is the square root branch point. Using equations (6.6) and 


(6.9) with the numerical values of a_p. 2 j, j = 0,1,2,..., obtained as described in the 
beginning of this Section, one finds in the first sheet a numerical value of 2:i,+,s«m(0) for 
each numerical value of Xc- Then numerical Newton (secant) iterations are performed 
over Xc aiming to ensure that zi^+^sum{0) converges to ~ Zpade{0), i-O. Xc is chosen such 
that Zi^+^sumiO) recovers the value of Zpade(0). It provides Xc with the precision at least 
10“^® which is limited by the precision of Fade approximation. Part I also demonstrated 
the calculation of Stokes wave well beyond quad precision by using variable precision 
arithmetics with the achieved accuracy ~ 200 digits thus increasing of accuracy for Xc 
is also possible if needed. Table of Appendix [C] provides numerical values of Xc which 
correspond to Fade approximations of Stokes wave found in Part I. 


6.2. Finding o/a+^ 2 n-i-i,o o,nd a_, 2 n,o from matching the series expansions at C = ±iXc 

in the second, third etc. sheets 


The procedure for finding numerical values of Xi 
Section 


6.1 


and a+pj-, j = 0,1, 2,... described in 


together with equations (6.3) and (6.7), allows to immediately find a_pj, 


j = 1,2,... for each given value of a-^ 2 ,o- Similar to equations (6.8) and (6.9), a notation 
is used such that ^:j,+,siim(C) ^;,-,siim(C) oxe the finite sums corresponding to the 

truncation of the series 2i,+ (C) (6.1) and Zi,-(C) ( |6.2[ ), respectively. We assume that 
jmax — 200 for these finite sums in the sheets I = 1,2,.... The numerical Newton 
iterations at the first step are performed over a-^ 2,0 aiming to ensure that Z 2 ,-,sum( 0 ) 
converges to Z2^+^sum(0). At the second step, the Newton iterations allow to find a+p_o by 
matching Z 3 _ sum(0) and 23 ,_,sum( 0 ). In a similar way, the third, fourth etc. steps allow 
to find a_p_ 0 ) a-i-, 5 , 0 ) a-.e.Oi Then using equations (6.7), one obtains values 


of a+^„_o and a-^n,o for all positive integer n completing the analytical continuation of 
Stokes wave into the disks ± iyd < 2xc in the infinite number of sheets of Riemann 
surface. 

The result of that analytical continuation was compared with the analytical continua¬ 
tion by ODE integration of Section [^giving the excellent agreement which is only limited 
by the standard numerical accuracy ~ 10“^® of Stokes wave in the physical sheet. In¬ 
creasing that accuracy of analytical continuation is straightforward by increasing jmax 
for the finite sums zi^+^sumiC) and z;,-,sMm(C) (and similar by increasing the accuracy 
for ODE integration) provided Stokes wave precision is increased. Table of Appendix 
[C] provides a sample of numerical values of a-^ 2 n,o, <x+^ 2 n+i,o, for n = 1, 2,3 obtained by 
the numerical method outlined in this Section. 


7. Singularities of Stokes wave for finite values of w 


Grant (1973) and Tanveer (1991) showed that the only possible singularity in the finite 


complex upper half-plane of the physical sheet of Riemann surface is of square root type. 
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This result is consistent both with the simulations of Part I ( Dyachenko et al^|2016 ) and 
numerical integration of ODE (3.4) in Section]^ 

The analysis of jTanve'CT (1991) is based on a version of equation (3.2) together with 
the assumption of the analyticity of z{w) in C'*' for the fist sheet of Riemann surface. 
Assume that one performs ODE integration in the second, third etc. sheets of Riemann 
surface as described in Sections and with z{w) at the nth sheet coupled to z{w) 
in the n — 1th sheet. Here the counting of sheets follows Section [5] and ass umes that 
—IT < Re{w) < TT, \Im{w)\ < oo for all sheets. Then the analysis of Tanveer ( 1991[ ) can 
be immediately generalized to the nth sheet at values of w = wi such that z{w) has 
no singularity a,t w = wi in the n — 1th sheet (see equation ( |7.16[ ) below). Coupling of 
square root singularities at C = ±ixc which is studied in Section [6| however goes beyond 
the analysis of Tanveer (|1991 ). 


The series expansions of Section show that square root singularities can occur at 
any finite values of w = Wi away from the real axis. It was found in Section that each 
square root singularity can either have a sister square root singularity at the complex 
conjugated point w = wi in the same sheet or can exists without the sister singularity at 


w = wi thus going beyond the case analyzed by Tanveer (1991). A question still remains 


if any other type (beyond square root) of coupled singularities at w = wi and w = wi in 
the same sheet is possible. 

Going from the first sheet to the second one, then from the second one to the third 
one etc., one concludes that the only way for the singularity other than square root to 
appear is to be coupled with the square root singularity in the previous sheet. Otherwise 
if would violate the above mentioned generalization of the result of [Tanve'er (1991) to 
the arbitrary sheet. Consider a general power law singularity of z(w) at w = wi coupled 
with the square root singularity of z(w) at w = wi. We write that general singularity in 
terms of double series as follows 


j(w) = 




(7.1) 


where a is the real constant, Cn,m are complex constants and n, m are integers. By shifting 
n and m one concludes that without loss of generality one can assume that 


0 < a < 1/2. 


(7.2) 


After the complex conjugation, the square root singularity of z(w) at w = Wi is given by 
the following series 


^(ra) = ^ dn(w - 


Wi 


,n/2 


(7.3) 


n=0 


where are the complex constants. Coupling of z(w) and z(w) in equation (3.3) explains 
why half-integer powers n/2 must be taken into account in equation ([TT])- It is convenient 
to transform from w into a new complex variable 


q=(w-Wi)^/^. 


Equation (3.3) for the new variable q takes the following form 

4g^c^ 


Zo = 


Zq[l{z - Z) + C^y 


(7.4) 


(7.5) 


The series (7.1) is then transformed into 
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,n+2ma 


(7.6) 


while the series (7.3) runs over integer powers, 


Kq) = '^dnq'^. 

n—0 


(7.7) 


The series (7.6) can be also called by -^-series, see e.g. Hille (1997). If a is the rational 


q thus reducing equation (7.1) to Puiseux series 


number then in equation (7.1) one can gather together all terms with the same power of 

(7.8) 


'(i) = X! 

n——(x> 


2njk 


where k is the positive integer. 

If one additionally restricts that there is no essential singularity at g = 0 then one has 


to replace (7.6) with the truncated series 


ziq) = (7.9) 

n^TiQ ,m^mQ 

for the integer constants ng and mg. Plugging in equations dTTl ) and ( |7.9[) into the Stokes 
wave equation (7.5), moving the denominator to the l.h.s. in equation (7.5) and collecting 


terms with the same power of g, starting from the lowest power, one obtains that 2a must 
be integer for any values of ng and mg and all values of dn- Thus no new solutions in the 
form (7.9) exists beyond what was found in Section]^ 

One can also study singularities using the classification of movable and fixed singular¬ 
ities in nonlinear ODEs of 1st order in the general form Zq = /(g, z) (see e.g Golubev 


(1950); Hille (1997); Ince ( 1956| )). Position of fixed singularities for the independent com¬ 
plex variable g is determined by the properties of ODE, i.e. by singularities of the function 
/(g, z). In contrast, the position of movable singularity is not fixed but typically is de¬ 
termined by an arbitrary complex constant. To analyze singularities, it is convenient to 
introduce a new unknown 


m ^ T 


\{z - z) + c 


,2 • 


Then equation (7.5) takes the following form 


(7.10) 




Aiq^c^Sf’ 


(7.11) 


where one reminds that z(g) is assumed to be known and is determined by z(g) from 
the previous sheet of Riemann surface. Equation (7.11) has a cubic polynomial r.h.s in f 
which ensures that it has a movable square root singularity 


OO 

e ^ c„(g - C)”/^ (7.12) 

n——l 


provided C 7 ^ 0, Zq[C) 0 (see e.g Golubev 

(1950) 

Hille 

(1997 

I; 

Ince 

(1956 

)), where 

Cn and C are the complex constants. Using equations ( 

7.4 

I, ( 

7.10 

1 and the condition 
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C ^ 0, one recovers the expansion (7.3) with replaced by wi + thus the movable 


singularity (7.12) is reduced to the square root singularity in w. 


Equation (7.111 has a fixed singularity at g = 0 provided Zq{0) ^ 0. To show that one 
uses a new unknown i/) = 1/^ to transform equation (7.11) into 


'^q = 


—iZqijj + 


Zqi’ 


(7.13) 


which has 0/0 singularity in r.h.s. for g = "0 = 0 satisfying the criteria for the existence 


of fixed singularity (see Golubev (1950); Hille (1997)). 

Consider now a particular case Zq{0) = 0 and Zqq{0) ^ 0 which corresponds to the 


expansion (6.3). One can define a new function 


5(9) ^ 


9(0) ^ 0 


(7.14) 


and rewrite equation (7.13) as follows 


i’q = 


—igg(g)^V' + 4igc^ 
9(9)V’ 


(7.15) 


Generally this equation still has a fixed singularity because of 0/0 singularity in r.h.s.. 
However, in a particular case when g{q) is the even function of g, i.e. one can define the 
function g(g^) = g{q) which is analytic in the variable g = g^ at g = 0. This case means 


that z{w) is analytic at w = wi. Then one transforms equation (7.15) into the equation 


tpq = 


—\g{qY'ij) + 4ic^ 

29(9)V' 


(7.16) 


which does not have a fixed singularity. Equation (7.16) together with equation (7.12) 


reproduces the result of Tanveer (1991) applied to all sheets of Riemann surface. 


Thus the approaches reviewed in Refs. Golubev (1950); Hille (1997); Ince (1956) ap¬ 


plied to equation (7.5) are consistent with square root singularities and the series expan¬ 


sions of Section for all sheets of Riemann surface. However, these approaches cannot 
exclude the possibility of existence of other types of singularity. Note that examples given 
Golubev (1950); Hille (1997|; Ince (1956) also show that the existence of the fixed sin¬ 


gularity in ODE at the point g = 0 does not necessary mean that the singularity occurs 
in the general ODE solution z{q) at that point. 

One concludes that a coupling of the essential singularity at w = Wi with the square 
root singularity at w = wi cannot be excluded neither by the series analysis used in 


equations (7.4)-(7.9) nor by looking at the fixed ODE singularities through equations 
(7.10)-(7.15). However, the simulations of Section and series expansions of Section 
1^ clearly indicates the absence of any singularities beyond square roots in all sheets 
of Riemann surface for non-limiting Stokes wave in C € C\((—ic»,—i] U ([i,ioo)) (i.e. 
everywhere in the complex plane C except the branch cuts (—icx), —i] and [i,ioo)). In 
the first sheet the branch cuts (—icx),—i] and [i,ioo) are not significant as explained in 
the Introduction, and the only non-square root singularity exists at C = b see Section 
1^ It is conjectured here that non-square root singularities do not appear in all sheets of 
Riemann surface for ( S C\((—ioo,—i] U ([i,ioo)). 

As discussed at the end of Section singularities are possible at the boundary of the 
strip Re{w) = ±7r which corresponds to the branch cuts [i, ioo) and [—ioo, —i) in Q plane. 
However, these branch cuts are separated by the distance tt from the origin in w plane 
(or by the distance 1 in C plane) and they cannot explain the formation of the limiting 
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Stokes wave as Vc —t 0. The same is true for the singularity at re —>■ ioo (C —?► i) analyzed 
in Section m 


8. Conjecture on recovering of 2/3 power law of limiting Stokes wave 
from infinite number of nested square root singularities of 
non-limiting Stokes wave as Xc —t 0 

One concludes from Sections]^ andthat the only possibility for the formation of 2/3 
power law singularity (1.1) of the limiting Stokes wave is through the merging together 


the infinite number of square root singularities from different sheets of Riemann surface 
in the limit Uc —>■ 0. The total number of square root singularities could be either finite or 
infinite for Vc > 0, both cases are compatible with the expansions of Section]^ (although 
the infinite number of singularities appears to hold for the generic values of the expan¬ 
sion coefficients of Section . Both numerical ODE integration of Section and series 
expansions of Section reveal that the number of sheets of Riemann surface related to 
singularities at /' = ±ixc exceeds several hundreds for a wide range of numerical values 
10“^ ^ Xc ^ 0.2. It suggests that the number of sheets is infinite for all values of Vc- In 
any case, the number of singularities must be infinite as Vc —t 0 . 

Here a conjectured is made that the limiting Stokes wave occurs as the limit Uc —> 0 
of the following leading order solution 


2 = i Y + cixy®\/C-iXc 




(C - + (-2ixc)^^^ \laiXc^ + \/iC - + (-2iXc)^/^ 


1/4 


\ 


+ \l + \\ otiXc'" + x! i.C - iXc)^/^ + (- 2 %)^/^ 


1/8 


1/4 


\ 


1/22„+2 

Cl2n+lXc 


\ 


+ \J ■■■ + ]/ «iXc^'‘ + \/(C - iXc)^/^ + (-2ixc)^/^ 

+ (C-iXc)^/^ + (-2iXc)^^^ ]J ocxxJ^ + -/(C - iXc)^/^ + (-2ixc)^^ 



X • • • -l-h.o.t., 


■ • ■ + Y “1^=^ + V “ iXc)^/^ + (-2ixc)^/^ 

( 8 . 1 ) 


which is the infinite product of increasingly nested square roots. This conjecture was first 


presented at the conference talk Lushnikov et al. (2015). Equation (8.1) has two terms 


with nested roots, one with nonzero complex constants ai, a 2 , cts, 04 ,... and another 
with nonzero complex constants di, 0 : 2 , a^, 0 : 4 ,... which are related through complex 
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conjugation as follows 


ai = a\e a-i = oie ..., a„ = a^e , 


( 8 . 2 ) 


All these constants including another complex constant c\ are of the order 0(1) inde¬ 


pendent of Xc- The relations (8.2) ensures that the symmetry condition (5.2) is satisfied. 

At C ^ Ac one obtains from equation (8.1) using the asymptotic of products of all 
square roots that 

Z (X ^1/2-H1/8+1/32-H/128+... ^ ^2/3 j-g 3 ^ 


exactly reproducing the Stokes solution (1.1) while the term c\xJ^— ixc vanishes as 
Xc —>■ 0. For small but finite Xc the limiting Stokes solution O) is valid for Xc C 
as seen from equation (8.1). For C ~ Ij the higher order terms denoted by h.o.t. both in 
equations O) and (|8.1[), becomes important such as the term with the irrational power 


oc C^, Ai = 1.4693457... 


(8.4) 


( lGrant|[T973l |Williams|[IM| ) . 


Different branches of all nested square roots in equation (8.1) choose different sheets 


of Riemann surface following the numeration of sheets used in Section In particular, 
the principal branch of (C — iXc)^^^ in the expression = (C, — iXc)"^ + (—2ixc)^^^ 
corresponds to the first sheet. To understand that one expands ( 7 (C) at C = ~iXc which 
results in 

C + iXc 


ff-h(C) = 2(-2ixc)^/^ + 


+ 0(C + iXc)^ 


(8.5) 


2(-2iXc)i/^ 

where the subscript “-I-” means taking the principle branch of (C~iXc)^^^- For the second 
(negative) branch of (( — ixc)^^^ one obtains that 


9-(0 = - 


C + iXc 


2(-2iXc)F2 


0(C + iXc)' 


( 8 . 6 ) 


with the subscript “ — ” meaning that second branch. 

The expression g(() enters under the most inner square root into each term of the 
product in equation EL Then using the expansion ( |8.5| ) one obtains that the series 
expansion of equation (8.1) at C = —iXc contains only nonnegative integer powers of 


C+iXc thus confirming that 2(0 is analytic at C = —iXc- It means that the condition (6.5) 


is satisfied. In contrast, taking the expansion (8.6) one obtains that the series expansion 


of equation (8.1) at ( = —iXc is of the type (6.2) containing nonnegative half-integer 
powers of C + iXc as expected for all sheets starting from the second sheet. In addition. 


the term g(() in the square brackets in equation (8.1) ensures that a__ 2 n.i = 0 as required 
by equation (6.3). The expansion of equation (8.1) at ( = ixc bas half-integer powers of 
C — ixc for both branches of g(() thus being consistent with the square root singularity 
at C = iXc in all sheets of Riemann surface including the first sheet in agreement with 
equations (|6.1[), (6.4) and (|6.6|). 


Choosing two possible branches of all other nested square roots (besides the most inner 
square root) in equation (8.1) one obtains the expansions (6.1) and (6.2), at C = ±iXc with 
different coefhcients and at each Ith sheet. The values of these coefficients are 
determined by both values of Xc Ci, Oi, 02 , <^ 3 , 04 ,. •. together with the contribution 


from h.o.t terms in equation (8.1). One concludes that the ansatz (8.1) is consistent with 


the properties of non-limiting Stokes wave studied in this paper which motivates the 


conjecture (8.1). Also coefhcients ai, 02 , 0 . 3 ,... determine additional square root branch 


points which are located away from the imaginary axis at distance larger than Xc from 
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the origin in the third and higher sheets of Riemann surface. Values of these coefficients 
can be determined from the locations of branch points thus independently recovered from 
ODE integration similar to described in Section]^ (see also Section 8.1 for the example 
of recovering of ai). 


8.1. Numerical verification of the conjecture 


Now we provide a numerical demonstration of the efficiency of the conjecture (8.1 1 by 


using the simplest nontrivial approximation of equation ( 8 . 1 ) which takes into account 
only the three-fold nested roots as follows 


+ ciXc^®\/C-iXc + caC 


C3X 


1/24 


(C - iXc)^/^ + (-2iXc)^''^ y+ \/(C - iXc)^/^ + (-2ixc)^ 


cae 8 


/1/24 


(C-iXc)'/" + (-2iXc)'/^ 


1/4 


aie 4 xc' ^ + (-2ixc)^/^ 

(8.7) 

where we added the term c^f as well the constants cq, C2 and C3 to approximate the 
neglecting of other nested roots (which includes 0 : 2 , as,...) in comparison with equation 
(8.1). In other words, we approximate all more than three-fold nested roots in equation 


(8.1) by Taylor series expansion keeping only the constant and linear terms in f in that 


expansion. The constant terms result in adding the constants C 3 and C 3 in front of nested 
roots as well as in the small correction cqXc^^ to the constant term ic ^/2 of equation 


(8.1)). The linear terms result in the appearance of the constant ci replacing higher 


order nested roots of equation |8.1[ ). Also the factor 036 s’" ensures the symmetry (5.2) 
and the additional scaling provides for C Xc the same total scaling oc x^J^ of the 


nested square roots of equation (8.7) as in equation (8.1). The approximation (8.1) can 


be valid only up to moderately large values of Re{C,) = f {a. comparison with equation 

( 8 . 1 ) suggests that it could be valid up to values of ( in about several tens of Xc after 


that higher order nested roots must come into play). 

To find the numerical values of ai, cq, ci, C2 and C3 without any fit we use the following 
procedure. At the first step we determine from the contour integration of Section the 
location (5.3) of the first pair of off-axis square root singularities (located at C = Cci 
and ( = —Cci)- Then we solve equation aix^^ + \/((a — iXc)^/^ + (—2ixc)^/^ = 0 
(corresponds to the zero under three-fold square root in equation (EH)) in the third 
sheet of the Riemann surface which together with (5.3) gives that 

ai ~-0.0955383 -i 1.8351. ( 8 . 8 ) 

Note that the second square root = —Cci (symmetric with respect to the imaginary 
axis) is ensured by the similar term with in equation (8.7). At the second step we 


expand equation (8.1) in the first sheet of Riemann surface in powers of (C~iXc)^^^- After 
that we match the first five coefficients of that expansion to the analytical expressions 
of the coefficients j = 1,..., 5, of the expansion ( 6 . 6 ) obtained in Section]^ That 

matching results in the explicit expressions for cq, Ci, C 2 and C 3 . E.g., for Stokes wave 
with Xc = 2.9691220994... • 10“^ (corresponds to the last line of table[^of Appendix 
see that Appendix for more details on the numerical Stokes wave) we obtain that cq = 
i 17.1920 ..., Cl = e-^’"/‘^3.81499..., C 2 = -1.8779... and C 3 = 1.42696 ... - i 1.8849... 
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Figure 4. A comparison of the limiting Stokes wave (dotted line) with equation ( |8.7[ ) and 
the numerical solution for Stokes wave (last two of these are shown by the single solid line 
because they are visually indistinguishable with maximum difference between them ~ 4 • 10~®) 
for Xc = 2.9691220994... • 10“^. Solid line corresponds to —50xc ^ C ^ 50xc- 


Changing of Xc by several orders of magnitude results in changing these coefficients only 
within the range 5 — 10%. 

We demonstrated the efficiency of the obtained numerical values of Cq, Ci, C 2 and C 3 
in two independent ways. In the first way, it was checked that the coefficients a^ i j of 


equation ( | 6 . 6 | ) for j = 6 ,... are well reproduced (within 4% and 7% accuracy for j ^ 10 
and j ^ 100, respectively) by the expansion of equation (8.7) with the same numerical 
values of cq, ci, C2 and C3. It implies that the approximate expression (8.7) captures 
the significant property of the convergence of the series ( 6 . 6 ) rather than being just the 
match of a few first terms of that series. The second way of efficiency demonstration 
is provided in Figure where the excellent agreement is shown between the numerical 
solution of Stokes wave from Part I ( Dyachenko et aL|2016 ) and the expression (8.7) for 
— 50xc ^ C ^ 50xc- That range of C is far beyond the disk of convergence |C ~ iXc| < “^Xc 
of the series ( 6 . 6 ). 


9. Concluding remarks 

In summary, it was found that the Riemann surface corresponding to non-limiting 
Stokes wave consists of the infinite number of sheets corresponding to the infinite number 
of square root branch points located at ^ = ±ixc in all sheets except the first sheet. The 
first (physical) sheet has only one singularity at C = iXc while avoiding singularity at 
^ = —iXc which ensures that Stokes wave represents the analytical solution inside the 
fluid domain. Two ways of analytical continuation into all these sheets were used, the 
first one is based on ODE integration of Section and the second one is based on the 
coupled series expansions ( | 6 . 1 [ )-( | 6 ?^ in half-integer powers at ^ = ±ixc- 

To go beyond the disks of convergence |C ± iXcl < 2xc of the series expansions (6.1 )- 
( 6 . 6 ), it is conjectured in Sectionj^that the leading order form of non-limiting Stokes wave 
consists of the infinite number of nested square roots (8.1). These nested square roots can 
recover the series expansions ( | 6 . 1 [ )-( [ 6 ^ within their disks of convergence |C±iXc| < 2 xc- 
For 1^ ± ixd ^ Xc well beyond these disks of convergence, the asymptotic (8.3) is valid 
thus ensuring that the nested square roots form 2/3 power law singularity of the limiting 
Stokes wave in the limit Xc 0. 

There are two other infinite sequences of Riemann sheets resulting from (a) off-axis 
square root singularities in the third and higher sheets of Riemann surface as analyzed in 
Sections]^ andand (b) the singularity at C = i (corresponds to w = ioo) which involves 
logarithms as analyzed in Section However, these extra sheets do not contribute to 
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the qualitative change of power law singularity from 1/2 (non-limiting Stokes wave) to 
2/3 (limiting Stokes wave) near the origin as given by the asymptotic (8.3). However, 
these extra sheets are expected to be important for the analysis of Stokes wave for 
~ 1, where the higher order terms becomes important such as the term (8.4) with the 
irrational power ( Grant|[l973 William^|1981 1. The analysis of these terms is beyond the 
scope of this paper. These terms might be also essential to answer the question left open 
by Refs. Longuet-Higgins & Fox (1977); McLeod (1997) about whether the number of 
oscillations in the slope of a non-limiting Stokes wave wave increases to infinity as non¬ 
limiting Stokes wave approaches its limiting form. Note that these oscillations vanish for 
the limiting Stokes wave as was proven in Ref. Plotnikov & Toland (2004). 


Appendix A. Equivalence of two form of equation for Stokes wave 


In this Appendix we show that both forms (2.1) and (3.1) of equation for Stokes wave 
are equivalent to each other. Then Section ^implies that equations (2.13) and (2.22) are 
also equivalent to equations (2.1) and (3.1). Equation (2.1) was obtained by Dyachenko 


et al. (19961 while equation (3.1) in slightly different forms was used by numerous authors 


including Grant (1973); Longuet-Higgins & Fox (1977); Schwartz (1974); Stokes (18806); 
Williams (1981) and Tanveer (1991). Appendix |B| explains the derivation of equation 
(3.1) starting from basic equations of the potential flow of ideal fluid with free surface. 


Applying the Hilbert operator H (2.2) to equation (2.1) and using the relations (2.14) 
and (2.15|) one obtains that 


c^x,. 


- yxu + II[yyu] = 0 , 


which is equivalent to equation (2.13). We define a new variable 


/ = -H[yyu] 


(Al) 


(A 2) 


and split it into two functions 

/ = /+ + /■, (A3) 

using equations (2.4), (2.7), (2.8) and (2.9) such that /+ and f~ are the functions which 
are analytic in upper half-plane and lower complex half-plane C“ of w, respectively. 


Note that the zeroth harmonic /o = 0 as it follows from the definition (A2|. Taking the 


linear combination of / and Hf, using equations (2.1), (Al) and (A2), one finds that 

XuHf Pynf = c^XuVu- (A 4) 

Using obvious relations 


1 / - ^ 1 . 

Xu — 2 T ) ') Uu — 2 ^ \^u ) 


together with equations (A3), (2.10) one obtains from equation (A4| that 


Zuf^ - Zuf = 



(A 5) 

(A 6) 


Recalling that is analytic in C and, respectively, Zu 
the projector (2.11) to equation (A6) and find that 


is analytic and in C"*", we apply 


/+ = 

r = 




4 Zu 


(A 7) 
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where 2 ^ ^ 0 in C and 0 in C"*" because z(w) is the conformal transformation in 
C-. 


Then using equations (All, (A2), (A3), (A5| and (A71 with some algebra we recover 
equation (3.1) thus completing the proof of its equivalence to equation (2.1). Note that 
the mean-zero elevation condition (2.151 is essential in that proof making equation (All 
valid. Shifting of the origin in y-direction would result in the nonzero value of the mean 

7T 

elevation ymean = ^ J r]ix,t)dx. Then one would have to replace y hy y — ymean in 


— n 


equation (3.1 

). E.g. 

Tanveer 

(1991 

1 took ymean = -C^/2. A similar choice ymean = -€^12 

was used by 

Grant 

1973 

1 ,^ 

A^illiams 

(19811 and 

Plotnikov ( 

1982 

) up to trivial scaling of 


parameters. 


Appendix B. Stokes wave in the rest frame and in the moving frame 

Starting from Stokes (1880&I, it has been common to write Stokes wave equation at 


moving reference frame in transformed form with a velocity potential and a stream func¬ 
tion used as independent variables, see e.g. 


Grant 

(1973) 

Williams 

(1981 

) and 

Tanveer 


(1991). The purpose of this Appendix is to relate that traditional form of Stokes wave 


equation to another form used by to Dyachenko et al. (19961; Zakharov et al. (20021. 


In physical coordinates (x, y) a velocity v of two dimensional potential flow of inviscid 
incompressible fluid is determined by a velocity potential ^{x,y,t) as v = Vd). Here x 
is the horizontal axis and y is the vertical axis pointing upwards. The incompressibility 
condition V • v = 0 results in the Laplace equation 


= 0 


(Bl) 


inside fluid —oo < y < r]{x,t). The Laplace equation is supplemented by the dynamic 
boundary condition (the Bernoulli equation at the free surface y = r]{x,t)) 


1h 


(V$)' 


-I- 77 = 0 


y=rj{x,t) 


and the kinematic boundary condition 


dr] 


dr] 9$ 
dx dx 


9 $ 

dy 


v='n(x,t) 


(B2) 


(B3) 


at the free surface. In our scaled units, the acceleration due to gravity is g = 1. We define 
the boundary value of the velocity potential as = ip(x,t). Equations 

(Bl), (B2) and (B3|, together with the decaying boundary condition at large depth 


^{x,y,t)\y^_ac = 0 


(B4) 


form the closed set of equations. Equation (|B 4|) implies that the rest frame is used such 


that there is no average fluid flow deep inside fluid. See also Part I (Dyachenko et al. 


2016) for more details on basic equations of free surface hydrodynamics. 


Consider the stationary waves moving in the positive x direction (to the right) with 
the constant velocity c so that 


$ = $(a; - cf,y), 
77 = r](x — ct). 


(B5) 


It was obtained in Ref. Dyachenko et al\ (1996) (see also Part I (Dyachenko et al. 
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20161 ) that 'h = —cHy = cx, where H is the Hilbert transform (2.2). Respectively, 
H'^ = cy. The complex velocity potential H at the free surface is given by 


n = 'I' + = c{x + yy — u). 


(B6) 


The analytical continuation of (B 6) into the lower complex half-plane w £ C is given 

by 


n = c{z — w) = cz. 


(B7) 


We perform a Galilean transformation to a frame moving with the velocity c in the 
positive direction with the new horizontal coordinate x' = x — ct so that the velocity 
potential and the surface elevation turn time-independent as $(a;') and rj{x'), respectively. 
Alternatively, one can also dehne a velocity potential in the moving frame as ^{x' ,y) = 
$(a; — ct, y) such that 


d) = (x — ct)c + $(x — ct, y). 


Then equation (B 2) results in 


y=r](x-ct) 


v-j]=o 


and (B3) gives 


dr] 

dx dx 



= 0 . 


(B8) 


(B9) 


(B 10) 


y = y{x-ct) 


The decaying boundary condition (B 4) is replaced by 


Equations (B9), (BIO) and (Bll 


moving frame, see e.g. Grant (19731; 


l>(x', j/)|y_,._oo = -c. (B 11) 

are the standard equations for Stokes wave in the 


Williams (1981). Often small variations of equations 


(B9), (BIO) and (Bll) are used such as a trivial shift of the origin in the vertical 


direction rj — ^—) 
and rescaling c to 


(1973). 


ry, assuming that Stokes wave moves in negative direction (to the left) 
one (then the spatial period 27T is also rescaled) as was done in Grant 


Similar to (B8|, we define the stream function in two forms, 0(x') and 0(x') (in the 
rest frame and in the moving frame, respectively) as follows 


0 = cy -I- 0(x — ct, y). 


(B12) 


Using equations (B8) and (B 12) one obtains that correspondingly, that two forms of the 

(B13) 


complex velocity potential, n(x') and n(x'), are given by 


n = d* -I- i0 = cz — c t -I- di -I- 10. 


A comparison of (B 7) and (B 13) reveals that 


n = d> -I- 10 = —c{w — ct) = —c 


(B14) 


where w' = w — ct. Thus H is the same as w' (up to the multiplication on —c) which 
explains why using the velocity potential d> and the stream function 0 as independent 


variables in Refs. 

Grant 

(1973) 

Stokes 

(18806); William 

3(1981 

) is equivalent to using w' 

as the independent variable in Ref. 

Dyachenko et al. 

(] 

996 

). The difference between H 
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and n is reflected by the boundary conditions (B4) and (B 11) such that for B the fluid 


at infinite depth has a zero velocity while for H that velocity is —c in the a:-direction . 
A technical advantage of working with B instead of B in Ref. Dyachenko et al. (1996) 


is that the decaying boundary condition (B 4) allows to relate real and imaginary parts 


of B through the Bilbert transform for real values of w' as 0 = and 4> = —HQ. 


Equations (B9) and (B14| results in Stokes wave equation in the form (3.1) after we 


notice that 




Appendix C. Tables for numerical values of Xc for Stokes wave 

Table provides a sample of the dependence of the singularity position Xc on the 
scaled wave height H/X = H/{2'k) for Stokes wave. Numerical values of Xc are obtained 
by the numerical procedure described in Section |6.1[ The Bade approximants from Part 
1 |Dyachenko et aL||2016|) (these appro ximants are also available through the electronic 
attachment to Ref. Dyachenko et al. (2015a I and at the web link Dyachenko et al. 


(20156)) are used for each values of H/\. More values of Xc for different values of H/X 
are also available at the web link Dyachenko et al. ( 20156[ ). The accuracy of numerical 
values of Xc is at least 10“^® which is limited by the precision of Fade approximation. 

We chose parameters at 1st, 3rd and 5th lines of table to correspond Stokes waves 
with c = 1.03, 1.066 and 1.086, respectively (here the exact values of c are used). These 
three particular values of parameters correspond to three highest Stokes waves provided 
in table 1 of Ref. Tanveer ( 1991[ ). Tablej^reproduces these three highest waves from table 
1 of Tanveer (1991), where the position of square root branch point C, = iXc is recovered 
from the parameter (q of Ref. |Tanveer| ( |199T| ) as Xc = ~(1 + Co)/(l ~ Co)- Also H in Ref. 
Tanveer ( 1991| ) is the half-height of Stokes wave so it is divided by tt in table The 
comparison of tables and [^reveals that while all digits except the last one or two agree 
for two smaller values oi H/X, but the agreement looses one more digit with the increase 
of H/X. It is possible that Ref. Tanveer (1991) expected that loss of numerical precision 
because the number of digits provided in table 1 of Tanveer (1991) decreases with the 
increase of H/X. 

Table [^provides a sample of numerical values of a__ 2 n,o and a+_ 2 n-i-i,o for four different 
values of Xc corresponding to table These numerical values of Xc are obtained by the 
numerical procedure described in Section [6.2[ For brevity only 16 digits of the numerical 
precision are shown. 


REFERENCES 

Amick, C. J. & Fraenkel, L. E. 1987 On the behavior near the crest of waves of extreme 
form. Trans. Amer. Math. Soc. 299 , 273-298. 

Amick, C. J., Fraenkel, L. E. & Toland, J. F. 1982 On the Stokes conjecture for the wave 
of extreme form. Acta Math. 148 , 193-214. 

Babenko, K. I. 1987 Some remarks on the theory of surface waves of finite amplitude. Soviet 
Math. Doklady, 35 ( 3 ), 599-603. 

Buffoni, B., Dancer, E. N. & Toland, J. F. 2000 The sub-harmonic bifurcation of Stokes 
waves. Arch. Ration. Mech. Anal. 152 , 241-271. 

Buffoni, B. & Toland, J. F. 2001 Dual free boundaries for Stokes waves. C. R. Acad. Sci. 
Paris Sr. I Math. 332 , 73-78. 

Dyachenko, Alexander I., Kuznetsov, Evgenii A., Spector, Michael & Zakharov, 
Vladimir E. 1996 Analytical description of the free surface dynamics of an ideal fluid 
(canonical formalism and conformal mapping). Rhys. Lett. A 221 , 73-79. 

Dyachenko, Sergey A., Lushnikov, Pavel M. & Korotkevich, Alexander O. 2013a The 
complex singularity of a Stokes wave. JETP Letters 98 (11), 767-771. 























































Branch cuts of Stokes wave on deep water. Part II 


35 


Wave height H/X 

0.077390566513510100664367446945009 

0.10042675172528485854673515635249 

0.11396866940628458279840665192065 

0.12063157457100181211171486096916 

0.13046836752896146189584028585057 

0.13871124459012593791450261565795 

0.14003037735536232024327827857514 

0.14011096764402710691403135029555 

0.14015101306439164612988663680930 

0.14033404782061154512392085005894 

0.14051416938624427610421738297959 

0.14056584420653835444977911685203 

0.14070850110629620828789822957203 

0.14074703013044272779483720282718 

0.14075662532618050016439516401203 

0.14077748818517580368147808000934 

0.14080831525231916769272562321913 

0.14083140371280991872217783523764 

0.14085072731982411577531399667650 

0.14086825990337854565346642922133 

0.14087792765270709969236336933758 

0.14088586197110133631188309127224 

0.14089635109209977336909824577330 

0.14091001709910523062648751945506 

0.14091839307555128402812965695553 

0.14092032625051507844376744407087 

0.14092252442341776630057428860718 

0.14092514757875551525458131241416 

0.14092738180637770768107092780251 

0.14093056906823728426117769727974 

0.14093510137194143743061264048898 

0.14094119430696937198416665739014 

0.14094867821783188240349944668053 

0.14095352707479979419800954052129 

0.14095778935504595764411825281530 

0.14097009565718766875950104063752 

0.14098407663748727496462567878823 

0.14100153154854889551064171690484 

0.14103365111671204571809985597404 

0.14105431648358048728514606849313 

0.14105777885488320816492860225696 


Singularity position Xc 
0.22959283981280615879703284574991 
0.12126855832745608069685459720991 
0.071654598419719678169515049620847 
0.050466513002046555340085106251597 
0.022711769117183995733113183176661 
0.0030056373876010407473234354599642 
0.0007999065189780408034349632263817 
0.00069951386487208337732279647662665 
0.00065164210434348698577048482811606 
0.00045087566212961727243263909506818 
0.00028427822364922236690177980170163 
0.00024252541408812956956630147113284 
0.00014199627497457559254018017702833 
0.00011868402545440790157599298606945 
0.0001131402886276901411780810604808 
0.00010145173966680681771175565637662 
0.000085108686515454366575393860892637 
0.000073606496213860270898473095913984 
0.000064475962982549833303295412314089 
0.000056590609636696915098098733019878 
0.00005240769363924544328892679639685 
0.000049063815868419517646209932713057 
0.00004476805660136311962064510052487 
0.000039388011825703454833655362993012 
0.000036214071851881467799287017287358 
0.00003549506133290811208694741093295 
0.0000346837089035969690548283554112 
0.000033724196620161039218316518297236 
0.000032914454078339616366407901860458 
0.000031771329157192593064752326105744 
0.0000301703287220913256069400404687 
0.000028063945797678144251500481216356 
0.000025549865907771481807832323273915 
0.000023964796260036642282422099761643 
0.000022600407539173053002286018858435 
0.000018816656490602043348418618380363 
0.00001480968355336403686583695738714 
0.000010273655389226364040855903301072 
3.5012288974834512273437793255939e-6 
6.0520035443913536064479745209207e-7 
2.9691220994639291094028846634237e-7 


Table 1. A sample of numerical values of Xc vs. the scaled Stokes wave height H/X. 


H/X c X. 

0.07739055 1.0300 0.22958 
0.1139758 1.0660 0.071667 
0.13055 1.0860 0.022769 


Table 2. Parameters of three highest Stokes waves of Table 1 from Tanveer (19911. Units are 
converted to the notation of this paper with the same number of digits kept as in lief. |Tanveer| 
IIQQlIl 
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Xc = 0.12126... Xc = 0.05046 ... 

a_, 2 ,o 1.947517181530394 1.332875450393561 
a+,3,0 1.933089192507101 1.395722669719572 
a_, 4 .o 2.823744469669705 1.830765178354924 
a+,5,0 2.715883020102187 1.841263995239744 
a_,6,o 3.541346294820654 2.238582675537623 


Xc = 0.000242... Xc = 2.969 ... X lO”’’ 

0.616114648091185 0.5967616372529635 

0.6227276830074182 0.5968472222666076 
0.630550992725188 0.5969268580437934 

0.6356646908933044 0.5969952862738779 
0.642377214720984 0.5970622067844041 


Table 3. A sample of numerical values of a_,2n,o and a+,2n+i,o, n = 1,2, 3, for different Xc. 
More accurate numerical values of Xc can be recovered from table 
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